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Abstract: 

We extend the lowest-order matching of tree-level matrix elements with parton showers 
to give a complete description at the next higher perturbative accuracy in a s at both small 
and large jet resolutions, which has not been achieved so far. This requires the combination 
of the higher-order resummation of large Sudakov logarithms at small values of the jet reso- 
lution variable with the full next-to-leading order (NLO) matrix-element corrections at large 
values. As a by-product, this combination naturally leads to a smooth connection of the NLO 
calculations for different jet multiplicities. In this paper, we focus on the general construction 
of our method and discuss its application to e + e~ and pp collisions. We present first results of 
the implementation in the Geneva Monte Carlo framework. We employ iV-jettiness as the jet 
resolution variable, combining its next-to-next-to-leading logarithmic resummation with fully 
exclusive NLO matrix elements, and Pythia 8 as the backend for further parton showering 
and hadronization. For hadronic collisions, we take Drell-Yan production as an example to 
apply our construction. For e + e~ — > jets, taking a s (mz) = 0.1135 from fits to LEP thrust 
data, together with the Pythia 8 hadronization model, we obtain good agreement with LEP 
data for a variety of 2- jet observables. 
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1 Introduction 

Accurate and reliable theoretical predictions for measurements at collider experiments require 
the inclusion of QCD effects beyond the lowest perturbative accuracy in the strong coupling 
a s . This is especially important in the complex environment of the LHC, which requires 
precise predictions for a broad spectrum of observables. Higher-order corrections in a s are 
important to predict total cross sections and other inclusive observables. Exclusive jet observ- 
ables, such as jet-vetoed cross sections, require the all-order resummation of logarithmically 
enhanced contributions. For many observables, an accurate description across phase space 
demands a combination of both types of corrections. For experimental analyses to benefit 
from these advances, it is crucial to provide the best possible theoretical predictions in the 
context of fully exclusive Monte Carlo event generators. 

The goal of modern Monte Carlo programs is to provide a proper description of the physics 
at every jet resolution scale. This is the motivation for the by-now standard combination 
of matrix elements with parton showers (ME/PS). [1, 2] Here, the parton shower provides 
the correct lowest-order description at small jet resolution scales, where the resummation 
of large Sudakov logarithms is needed, while at large jet resolution scales the exact tree- 
level matrix elements are needed to provide the correct lowest-order description. Hence, the 
ME/PS merging provides theoretical predictions at the leading 0(1) accuracy relative to the 
lowest meaningful perturbative order. Once one has a consistent matching between these two 
limits of phase space, the possibility to include exact tree-level matrix elements for several 
jet multiplicities follows almost automatically by iteration. 

Given the necessity of higher-order perturbative corrections to make accurate predic- 
tions, it is important to extend the perturbative accuracy of the Monte Carlo description 
to a s accuracy relative to the lowest order. This requires including the next higher-order 
corrections that are relevant at each scale. At small scales, i.e., small values of the jet resolu- 
tion variable, this requires improving the leading-logarithm (LL) parton shower resummation 
with higher-order logarithmic resummation, while at large scales this requires including the 
fully differential next-to-leading order (NLO) matrix elements. It is important to realize that 
typically a large part of phase space, often including the experimentally relevant region, is 
characterized by intermediate scales, i.e., by a transition from small to large scales. In the 
end, providing an accurate description of this transition region requires a careful combination 
of both types of corrections. 

Such a Monte Carlo description at relative 0(a s ) accuracy across phase space has never 
been achieved and is the subject of our paper. (We briefly summarize the existing efforts 
to combine NLO corrections with parton showers in section 1.1 below.) The crucial starting 
point in our approach is that all perturbative inputs to the Monte Carlo are formulated in 
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terms of well-defined physical jet cross sections [3, 4]. This allows us to systematically increase 
the perturbative accuracy by incorporating results for the relevant ingredients to the desired 
order in fixed-order and resummed perturbation theory. 

An essential aspect of any higher-order prediction is a reliable estimate of its perturbative 
uncertainty due to neglected higher-order corrections. To the extent that parton shower 
Monte Carlos provide perturbative predictions, they should be held to the same standards. 
An important benefit in our approach is that we have explicit control of the perturbative 
uncertainties and are able to estimate reliable fixed-order and resummation uncertainties. As 
a result, in Geneva each event comes with an estimate of its perturbative uncertainty, i.e., 
Geneva provides event-by-event theory uncertainties. 

In our approach, the Monte Carlo not only benefits from the resummation, but in turn also 
provides important benefits to analytic resummed predictions. For one, it greatly facilitates 
the comparison with experimental data, as it allows easy application of arbitrary kinematic 
selection cuts, which can often be tedious to take into account in analytic predictions. More 
importantly, resummed predictions require nonperturbative inputs which can be treated as 
power corrections at intermediate scales but become 0(1) corrections at very small scales. 
Here, these are provided "on-the-fly" by the nonperturbative hadronization model. In essence, 
we are able to combine the precision and theoretical control offered by higher-order resummed 
predictions with the versatility and flexibility offered by fully exclusive Monte Carlos. 

In this paper, we focus on the theoretical construction. We leave a discussion of the 
implementation details of the Geneva Monte Carlo framework to a dedicated publication. 1 
We will however highlight some of the main technical issues we had to overcome and discuss 
some implementation details in the application sections. In the remainder of this section, 
we briefly summarize the existing efforts to include NLO corrections in parton shower Monte 
Carlos and give a short overview of our basic construction. In section 2, we discuss in 
detail the requirements to obtain full a s accuracy as well as our method to achieve it. In 
section 3, we discuss the application to e + e~ — > jets, where we combine next-to-next-to- 
leading logarithmic (NNLL) resummation with NLO matrix elements, and present results 
from the implementation in Geneva together with a comparison to LEP measurements. In 
section 4, we discuss the application to hadronic collisions and show first results for Drell-Yan 
production, pp — > Z/j* — > £ + l~+ jets, obtained with Geneva. We conclude in section 5. 

1.1 Previous Approaches Combining NLO Corrections with Parton Showers 

Over the past decade, many steps have been taken to include NLO corrections into Monte 
Carlo programs [5-21]. By now, the MC@NLO [9, 10] and Powheg [15-17] methods are 
routinely able to consistently combine the fixed NLO calculation of an inclusive jet cross 
section for a given jet multiplicity with additional parton showering. These methods have also 
been extended to include the full tree- level matrix elements for additional jet multiplicities [3, 
22, 23]. 

lr The current Geneva framework and implementation is new and independent of the earlier work in Refs. [3, 

4]- 
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Recently, efforts have been made to extend these approaches in order to combine NLO 
matrix elements for several jet multiplicities with parton showers [24-29]. We will discuss 
some issues faced by some of these approaches in section 2.1.5. Here, we would like to stress 
that including several NLO matrix elements by itself does not provide a full extension of 
the lowest-order ME/PS matching to relative 0(a s ) perturbative accuracy, since the fixed 
NLO corrections only suffice to increase the perturbative accuracy in the region of large jet 
resolution scales. To the same extent that the inclusion of the LL Sudakov factors in the 
ME/PS merging are needed to get meaningful results at intermediate and small jet scales, 
higher-order resummation is necessary to improve the perturbative accuracy in this region. 

In our approach, the full information from NLO matrix elements for several jet multi- 
plicities is automatically included as follows: For a given Born process with N partons, a 
small jet scale corresponds to the exclusive iV-jet region, and here the 7V-parton virtual NLO 
corrections are incorporated in conjunction with the higher-order resummation; in fact, they 
are a natural ingredient of it. On the other hand, a large jet scale corresponds to the inclu- 
sive (N + l)-jet region with additional hard emissions. Here, the (N + l)-parton virtual NLO 
corrections are included in the usual way by the fixed NLO calculation for N + 1 jets. 

1.2 Brief Overview of Our Construction 

The starting point of our approach is the separation of the inclusive iV-jet cross section into 
an exclusive iV-jet region and an inclusive (N + l)-jet region: 

where da/d&N+i(T) denotes the differential cross section for a given value of T. Here T = 
T(&n+i) is a suitable resolution variable, which is a function of 3>Ar+i, and do /d<&N+i{T) 
denotes the fully differential cross section for a given T ■ In ME/PS merging, this role is 
played by the variable that determines the merging scale. The parameter 7~ cut serves as an 
infrared cutoff for the calculation of do /d<&N+i(T) and ideally is taken as small as possible. 

In the iV-jet region at small T (both above and below T" €Ut ), logarithms of T become 
large and must be resummed to maintain consistent perturbative accuracy to some order in 
a s . On the other hand, in the (iV + l)-jet region at large T, a fixed-order expansion in a s 
will suffice. To consistently match the resummed and fixed-order calculations, we use the 
following prescription for the jet cross sections entering in eq. (1.1): 



do J_resum 



d$AT d<3?jv 



da FO _ cut , do 



(X )-^ — (T c 



d$iv d<l?Ar 



FOJ 



do do 



FO 



d$7v + i d® N+1 



do msum I do 



rcsmn 



d^ArdT / d$ N dT 



FO 



(1.2) 



The superscript "resum" indicates an analytically resummed calculation and "FO" indicates 
a fixed-order calculation or expansion. This construction properly reproduces the fixed-order 
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calculation at large 7~, the resummed calculation at small T, and smoothly interpolates 
between them. 

It is straightforward to extend our formulation to combine higher jet multiplicities at NLO 
with higher-order resummation, as we will show. This is done by replacing dcr FO /d^^+i m 
eq. (1.2) with an inclusive (N + l)-jet cross section separated into the exclusive (N + l)-jet 
and inclusive (JV + 2)-jet cross sections and iteratively applying eq. (1.2). 

The key ingredients in our approach are the higher-order resummation of the jet resolution 
variable, the fully differential fixed-order calculation, and the parton shower and hadroniza- 
tion. While each of these components are known, there is a sensitive interplay of constraints 
between them that must be satisfied to achieve a consistent combination. This is precisely 
what is accomplished in the Geneva framework and is the focus of this paper. 

2 General Construction 

In this section, we derive our theoretical construction in a process-independent manner. We 
start in section 2.1 with a slightly simplified setup, considering the singly differential spectrum 
in the jet resolution variable. We use this to discuss in detail the perturbative structure and 
the accuracy in the different phase space regions. In section 2.2, we discuss the extension 
to the fully differential case and how to combine the fixed-order expansion and resummation 
in this situation. In section 2.3, we further generalize these results to include several jet 
multiplicities by iteration. Finally, in section 2.4, we discuss the Monte Carlo implementation 
and how to attach parton showering and hadronization. 

2.1 What Resummation Can Do for Monte Carlo 
2.1.1 Basic Setup 

The basic idea of Monte Carlo integration is to randomly generate points in phase space 
("events") that are distributed according to some differential (probability) distribution. By 
summing over all points that satisfy certain selection criteria, we are able to perform arbitrary 
integrals of the distribution. In our case, that distribution is the fully differential cross section, 
allowing one to compute arbitrary observables. For simplicity, we will first focus on the singly 
differential cross section in some phase space resolution (or jet resolution) variable T of 
dimension one. The precise definition of T is not important at the moment, so we keep it 
generic for now. We use the convention that the limit T — > corresponds to Born kinematics, 
i.e., the tree- level cross section is ~ <5(T). We also require that T is an IR-safe observable, 
such that the differential cross section da/dT can in principle be well defined to all orders in 
perturbation theory and for T > contains no IR divergences. 

To give an example, for our application to e + e~ — > 2/3 jets in section 3, we will use 
2-jettiness 7i = E cm (l — T), where T is the usual thrust [30]. Alternatives include other 2-jet 
event shapes. For Drell-Yan in section 4, we will use beam thrust [31]. An alternative would 
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be the pr of the leading jet. If the Born cross section we are interested in has ./V signal jets, 2 
then T could be iV-jettiness or the largest pt of any additional jet. The important point 
is that we can think of T as a resolution variable which determines the scale of additional 
emissions in the <3?> n + i phase space, such that for T < 7" cut there are no emissions above the 
scale T cut . For later convenience, we also define the dimensionless equivalent of T as 



Here, Q is the relevant hard-interaction scale in the Born process, e.g., Q = E cm for e + e~ —> 
jets or Q = mp+g- for Drell-Yan pp — > Z/j* — > l + £~. In terms of r, the limit r <C 1 
corresponds to the exclusive limit close to Born kinematics. For t ~ 1, there are additional 
emissions at the hard scale T ~ Q, which means we are far away from Born kinematics and we 
should switch the description to consider the corresponding Born process with one additional 
hard jet. 

To describe the differential T spectrum, we want the Monte Carlo to generate events at 
specific values of T, which are distributed according to the differential cross section da/dT. 
The total cross section is then simply given by summing over all events, 



The essential problem every Monte Carlo generator faces is that in perturbation theory the 
differential cross section da/dT contains IR divergences from real emissions for T — > 0, 
which only cancel against the corresponding virtual IR divergences upon integration over the 
T — > region. As a result, the perturbative spectrum for T — > can only be defined in a 
distributional sense in terms of plus and delta distributions [see eq. (2.6) below]. To deal with 
this, we have to introduce a small cutoff 7~ cut and define the cumulant of the spectrum as 



In the Monte Carlo, the total cross section is then obtained by combining the cumulant and 
spectrum as 



In practice, this is implemented by generating two distinct types of events: (i) events that 
have T = and relative weights given by o"(T" €Ut ), and (ii) events that have nonzero values 
T > T cut and relative weights given by da / dT- The first type of events have Born kinematics 
and represents the tree-level and virtual corrections together with the corresponding real 
emissions integrated below T" cut . The second type of events contains one or more partons in 
the final state, since the real-emission corrections determine the shape of the spectrum for 
nonzero T. We now have two basic conditions: 

2 As usual, we assume that the Born cross section is defined with appropriate cuts on the iV signal jets, so 
that it does not contain any IR divergences by itself. 




(2.2) 




(2.3) 




(2.4) 
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1. From a numerical point of view, we want the value of T cut to be as small as possible, 
so as to describe as much differential information as possible. In practice, our ability to 
reliably compute the cumulant o(T cut ) in perturbation theory sets a lower limit on the 
possible value of T cut . 

2. Since T cut is an unphysical parameter, the dependence on it should drop out (to the 
order one is working). In practice, this is equivalent to including the corresponding 
dominant higher-order corrections in the cumulant and spectrum. 

2.1.2 Perturbative Expansion and Order Counting 

In perturbation theory, the differential cross section in r and the cumulant in r cut have the 
general form 



da 

d7 



dcr smg da non S 

+ 



dr 



dr 



a{r cut ) 



-f 

Jo 



dr 



da 
d7 



— (jSing^cut 



(r cut ) + a nons (T cut ) , 



(2.5) 



where we have distinguished "singular" and "nonsingular" contributions. For r — > 0, the 
singular terms in dcr smg /dr scale like 1/r, while the nonsingular terms in d<r nons /dr contain 
at most integrable singularities. For the cumulant, this means that cj smg (r cut ) contains all 
terms in cr(r cut ) enhanced by logarithms ln fc (r cut ), while ,j nons ( r cut = 0) = 0. 
The singular part of the spectrum is given by 



Cluing 

dr 



C^( as )5(T) + Y,Cn(a s )C n {r) 



(2.6) 



n>0 



where ob denotes the Born cross section, and we denote the usual plus distributions as 

~6{x)\n n {x) 



C n (x) 



I dxC n {x) 
Jo 



ln" +1 (x cut ) 



n + 1 



(2.7) 



They encode the cancellation between real and virtual IR divergences. The corresponding 
singular contribution to the cumulant cross section integrated up to r < r cut is 



sing/ cut 



(r cut ) = ob 



C_i(a s ) + ^2C n (a s 



n>0 



ln n+1 (r cut ) 
n+1 



(2.8) 



At O(cts), only the coefficients C n (a s ) with n < 2k — 1 contribute, so da/dr has logarithms 
up to a"L 2n_1 /r, while cr(r cut ) has logarithms up to a"L 2 " t , where we use the abbreviations 



L = lnfr) 



L cut = ln(r' 



cut\ 



(2.9) 



The a s expansion of the coefficients C_i(a s ) and C n {a s ) in the singular contributions 
can be written as 



C-i(a s ) = l + ^c fe _iaj, C n {a s ) = ^ c kn a[ 

k>l 2fc>n+l 



(2.10) 
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Figure 1. Illustration of the different parametric regions in the jet resolution. 



Similarly, the a s expansion of the nonsingular contributions can be written as 

dCT nons 

— = o B 

k>l 

Using eqs. (2.10) and (2.11), the spectrum and cumulant up to 0{a 2 s ) are given by 



°B £ fr S (r) a k s , Fr ns (r cut ) = f dr /^(r) . (2.11) 

J o 



1 da 
as dr 



T>0 



T 

2 



cnL + cio 



+ r/r ns (r) 



+ 



n 



T 



c 23 L 3 + c 22 L 2 + c 21 L + c 20 + r/ 2 nons (r) 



(2.12) 



0"B 



a ( T cut ) = l + a s — Leut + cio^cut + ci,-i 



+ «1 



C -fLi ut + + ^L c 2 ut + c 20 L cut + c,,^ + F 2 nons (0] + O(o») 



Note that the _i constant term in the singular corrections, which contains the finite virtual 
corrections to the Born process, only appears in the cumulant. 

We now distinguish three parametrically different regions in r, which are illustrated in 
figure 1: 

• Resummation ("peak") region t < 1: In this limit, the logarithms in the singular 
contributions are large, such that parametrically one has to count 3 



a s L ~ 1 , 



a sL cut ~ 1 . 



(2.13) 



In analytic resummation, the counting and resummation of logarithms is performed in the exponent of 
the cross section, where one counts a s L ~ 1. For the purpose of our argument in this section, it is sufficient 
to adopt the weaker scaling in eq. (2.13) and only count logarithms in the cross section. In our results, we 
always perform the full resummation in the exponent, as discussed in section 3.1.1. 
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This means one has to resum the towers of logarithms (a s L ) n in the spectrum and 
(a s L 2 ut ) n in the cumulant in eq. (2.12) to all orders in a s to obtain a meaningful per- 
turbative approximation at some order. At the same time, the nonsingular corrections 
can be regarded as power suppressed, since they are of relative 0(t). 

• Fixed-order ("tail") region r ~ 1: In this limit, the logarithms are not enhanced, and 
a fixed-order expansion in a s is applicable. The singular and nonsingular contributions 
are equally important and both must be included at the same order in ot s . In particular, 
there are typically large cancellations between these for r ~ 1, so it is actually crucial 
not to resum the singular contributions in this region, since otherwise this cancellation 
would be spoiled. 

• Transition region: The transition between the resummation and fixed-order regions. 

There are of course no strict boundaries between the different regions. This is why it is 
important to have a proper description not just in the two limits but also in the transition 
region, which connects the resummation and fixed-order regions. In fact, in practice the 
experimentally relevant region is often somewhere in the transition region, where both types 
of perturbative corrections can be important. 

2.1.3 Lowest Perturbative Accuracy 

For the Monte Carlo to provide a proper description at all values of T, it has to include at 
least the lowest-order terms relevant for each region. Keeping only these, and dropping all 
other terms, the spectrum and the cumulant are given by 



1 


da 




OB 


d7 


T>0 


1 


a(r cut ) 


OB 







T 



L h(a s L 2 ) + h(a s L 2 ) +r/ 1 nons (r) 



J cut/ 



(2.14) 

where the functions /o,i and i*b,i are given in terms of the coefficients Cij in eq. (2.10) as 

LL a : f (a s L 2 ) = £ c n+1 , 2n+1 (a s L 2 ) n , F (a s L 2 ) = ^ < 2 +1 f+ 1 (a s L 2 ut ) n , 
^-^ 2(n + 1) 

n>0 n>0 v ' 

NLL CT : h(a s L 2 ) = V c n+1 , 2n (a s L 2 )" , F^L 2 ) = V ^±^.(a 8 L 2 cut ) n . (2.15) 
^-^ ^-^ In + 1 

n>0 n>0 

The /o and Fq resum the leading-logarithmic series in the cross section, which we denote 
as LLq-. The functions fi and F\ resum the next-to-leading- logarithmic series in the cross 
section, which we denote as NLL CT . 

In the resummation region at r 1, the LL<j terms in the spectrum scale as L ~ l/y / o7 
(relative to the overall a s /r scaling) and provide the lowest level of approximation. The NLL CT 
terms scale as ~ 1, and one can argue about whether they are needed as well in order to get 
a meaningful lowest-order prediction. Formally, they are necessary to obtain the spectrum 
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at ~ a s / T , which one might consider the natural leading-order scaling of the spectrum (or 
equivalently if one does not want to rely on the ~ 1/ \fcts enhancement of the LL series). 
Experience shows that the NLL terms are indeed numerically important. For example, in 
analytic resummations, one rarely gets a sensible prediction without going at least to NLL. 
Similarly, to obtain sensible predictions from a parton shower, it is almost mandatory to 
include important physical effects such as momentum conservation in the parton splitting 
and the choice of a s scale [32]. In the cumulant, the LL^ series in Fq scales as ~ 1 and must 
be included. The NLL^ series in F\ scales as ~ y/a~ s and, for consistency, should be included 
in the cumulant if it is included in the spectrum. 

In the fixed-order region at r ~ 1, the lowest meaningful order in the spectrum is given 
by the complete 0(a s ) terms, requiring one to include the en and cio terms, which are part 
of the /o and f\ functions, as well as the nonsingular corrections /f ons (T). Since we take 
r cut to be small, the cumulant is always in the resummation region. Hence, its nonsingular 
corrections Ff° nB (r) [see eq. (2.12)] are suppressed by C(a s r cut ) and can be safely neglected. 

The leading level of accuracy in eq. (2.14) closely corresponds to what is achieved in 
the standard ME/PS matching. In this case, the LL resummation is provided by the parton 
shower Sudakov factors (either generated by the shower or multiplied by hand), where the 
jet resolution variable corresponds to the shower evolution variable, since that is the variable 
for which the shower directly resums the correct LL<j series. The LL CT series has a well known 
and simple exponential structure, 



Cn+l,2n+l 



r n+l 

2 n n! 



fo(a s L 2 ) = exp 



c ii T 2 
— ol s L 
2 



such that 



1 da 

a B dr 



T>0 



L 

en a s — exp 

T 



Cll 



1 



a(r cut ) = exp 



en 



2 a s-^cut 



(2.16) 



(2.17) 



The resummation exponent at LL CT is given by the integral over the leading en a s ln(r)/r term 
in the spectrum. This is precisely what the standard parton shower veto algorithm exploits 
to generate the resummation exponent. The analogous structure does not hold at NLL CT , 
which is why the parton shower cannot resum the NLL CT series by exponentiating the integral 
of the ciQa s /T term. As already mentioned, in practice, parton showers include important 
partial NLL effects, so practically this provides a numerically close approximation to the 
correct NLL CT series. The nonsingular corrections in the spectrum, /f ons (T), are obtained by 
including the full tree-level matrix element for one additional emission. Since the full matrix 
element also includes the en and cio terms, this requires a proper matching procedure to 
avoid double counting these terms. At LL CT , a simple way to do this is to multiply the full 
fixed-order result from the matrix element with the shower's LL CT resummation exponent, 



1 da 

as dr 



T>0 



T 



Cll L + c w + rfr S (r) 



exp 



c ii r 2 
— — a s L 

2 



(2.18) 
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which corresponds to the CKKW-L [1, 2, 33, 34] procedure. The reason this gives the spectrum 
correctly at LL<j is the simple structure in eq. (2.17), where the LL<j exponent multiplies the 
en term in the spectrum. 4 At large r ~ 1, the exponent in eq. (2.18) can be expanded as 
1 + (D(a s ), so eq. (2.18) gives the correct leading fixed-order result. 

Compared to eq. (2.14), the NLO matching performed in MC@NLO and Powheg 
amounts to adding to the cumulant the c\ _i singular constant, containing the 0(a s ) vir- 
tual corrections, as well as the nonsingular contributions F^ ons (r). Assuming the same set of 
NLL terms are included in the cumulant and spectrum, this achieves that inclusive quantities 
that are integrated over a large range of r, such as the total cross section, are correctly re- 
produced at fixed NLO, which provides them with 0(a s ) accuracy. In these approaches, the 
goal is not to improve the perturbative accuracy of the spectrum (or the cumulant at small 
r cut ), which has the same leading accuracy as in eq. (2.18). 

2.1.4 Next-To-Lowest Perturbative Accuracy 

We now want to improve the Monte Carlo description in eq. (2.4) from the lowest-order 
accuracy, given by eq. (2.14), to the next-to-lowest perturbative accuracy in a s . This requires 
us to include the appropriate higher-order corrections in each region, which gives 



1 da 
gb dr 



T>0 



r L 



+ 



a; 



-a it 



L/ (a s L 2 ) + / 1 (a s L 2 ) + r/ 1 nons (r) 

Lf 2 (a s L 2 ) + Ma s L 2 ) + rfr s (r 
1 + a s \L 2 cut F (a s L 2 cut ) + L cut iq(a s L 2 ut ) + ci,_i + F^ m {r mt ) 



+ a. 



L lut ^(as^cut) + Lcut F^{a s L 2 cnt ) 



(2.19) 



where we denote the series of logarithms resummed by the functions f 2 and F 2 by NLL^ and 
the series resummed by /g and Fg by NNLL<j. They can again be written in terms of the Cij 
coefficients in eq. (2.10) as 



nll; 

NNLL„- 



/ 2 (a s L 2 ) 
Ha s L 2 ) 



E 

n>0 

E 

n>0 



C n +2,2n+l 



(a s L 2 )\ F 2 {a s L 2 ) = Y, 



C n +2,2n+l , r2 s 



n>0 



2(n + l) 



C n +2,2n 



(a s L 



2\n 



F 3 (a s L 2 ) = 



n>0 



Cn+2,2n , t2 n 



(2.20) 



3/2 



and thus 



In the resummation region, the NLL^. series in the spectrum scales as ~ a 
provides the ~ a s correction to the LL CT series in /o. Similarly, the NNLLo- series scales as 
~ a 2 providing the ~ a s correction to the NLL CT series in f\. They can again be obtained 
by performing the standard resummation in the exponent of the cross section to NLL' and 



4 As before, since this simple LL CT structure does not hold at NLLo-, this procedure does not yield the 
resummed spectrum at NLLo-, even if one were to multiply the spectrum with the NLLo- resummation exponent 
as originally proposed in Ref. [1]. 
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NNLL respectively. (Here, NLL' refers to those parts of the full NNLL resummation that 
arise from the combination of the one-loop matching corrections with the NLL resummation, 
see section 3.1.1 and table 2.) 

In the fixed-order region, increasing the perturbative accuracy by ~ a s requires the 
complete O(af) corrections, including the f2° ns {T~) nonsingular corrections. Similarly, for the 
cumulant, F2 and F3 resum the NLL^ and NNLL CT series of logarithms, which scale as ~ a s 

3 II 

and ~ a s , respectively, and provide the ~ a s improvement over the LL CT and NLLo- series 
in Fq and F\. In addition, going to the next-higher order in the cumulant requires including 
the full singular constant c\—i as well as the nonsingular corrections Ff ons (T), which both 
scale as ~ a s . 

It is instructive to see where the information from the virtual NLO matrix elements enters 
in eq. (2.19). As already mentioned, the virtual NLO corrections to the Born process are given 
by c\-\. In addition, by multiplying the LL series it contributes part of fi and F2. Hence, 
consistently combining the virtual corrections with the resummation requires one to go to 
at least NLL'. The virtual NLO corrections with one extra emission (plus the integral over 
the two-emission tree-level matrix element) yield the full 0(a 2 s ) corrections in the spectrum, 
i.e., both the singular C2k terms as well as the nonsingular /2° ns terms in eq. (2.12). Adding 
these corrections again requires one to avoid double counting the singular C2k terms that are 
already included in the resummation. In analytic resummation, it is well known how to do 
this, namely by simply adding the nonsingular corrections. These are obtained by taking 
the difference of the full NLO corrections and the singular NLO corrections, where the latter 
are given by expanding the resummed result to fixed order. Since this construction involves 
the virtual contribution to both the Born process with one extra emission, we see that going 
consistently to higher order in both the resummation and fixed-order regions naturally leads 
to a combination of the information from two successive NLO matrix elements. 

2.1.5 Merging NLO Matrix Elements with Parton Shower Resummation Only 

We stress that for a description at the next-higher perturbative accuracy across the whole 
range in r, it is not sufficient to include the fixed NLO corrections to the spectrum and take 
care of the double counting with the parton shower resummation. This only provides the 
proper NLO description in the fixed-order region at large t. In the transition and resummation 
regions, a proper higher-order description necessitates higher-order resummation. Of course, 
this is not a problem if the only goal is to improve the fixed-order region at large r as is the 
case for example in a recent MC@NLO publication [28]. 

However, including the fixed NLO corrections outside the fixed-order region, as is done in 
Sherpa's recent NLO merging [26, 27], can actually make things worse in two respects: First, 
numerically this will typically force the spectrum to shift toward the fixed-order result and 
away from the resummed one. Since this can shift the spectrum in the wrong direction, it can 

5 Formally, 01,-1 belongs to the NLL4 series in the cumulant, but for the sake of discussion, we keep it 
explicit. 
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potentially make the result less accurate. 6 At the same time, the perturbative uncertainties 
from fixed-order scale variation decrease, which only aggravates this problem. Multiplying the 
NLO corrections to the spectrum with LL parton shower Sudakov factors (see, e.g., Ref. [35]) 
can mitigate this to some extent but does not solve the problem. The only consistent way 
to include the fixed NLO corrections to the spectrum outside the fixed-order region, and 
in particular obtain reliable perturbative uncertainties, is to properly combine them with a 
higher-order resummation. 

Second, this explicitly spoils the formal 0{a s ) accuracy of the inclusive cross section. 
To see this, consider adding the fixed NLO corrections to the lowest-order spectrum and 
cumulant in eq. (2.14), properly taking care of the double counting at 0(a 2 ), which gives 



1 do- 
er b dr 

1 



T 



L/ (a s L 2 ) + / 1 (a s L 2 )+r/ 1 nons (r) 



T>0 

a(r cut ) = 1 + a s 



+ 



T L 



c 2 iL + c 20 + Tfr S (r) 



Li ut F (a s Li ut ) + L cut F^asLiJ + Cl ,_ a + 



(2.21) 



Using these expressions yields for the inclusive cross section, 



-a 



-a(r cut )+ f dr 



OB 

1 + a s 



1 dcr 
cjb dr 
ns / 1 n 



ft* 



021 r2 , r 
"^"^cut + c 20-^cut 



(2.22) 



While the first two terms give the correct NLO inclusive cross section, the 0(a 2 ) terms 

3 /2 

induced by the fixed NLO corrections in the spectrum formally scale as a s and a s and 
therefore spoil the formal 0(a s ) perturbative accuracy for the inclusive cross section and in 
fact for any inclusive observable. This directly contradicts the claim in Refs. [26, 27] that this 
description maintains the higher-order accuracy of the underlying matrix elements in their 
respective phase space range. It only preserves the fixed 0{a s ) terms, which in the context 
of combining fixed-order corrections with a logarithmic resummation is necessary but not 
sufficient to preserve the higher perturbative accuracy. 

This problem cannot be avoided by multiplying the a 2 corrections in the spectrum with 
the LL parton shower Sudakov factors, since this does not provide the proper NLL' CT and 
NNLL CT series. Note also that we have already assumed in eq. (2.22) that the full NLL CT series 
is included in the spectrum and cumulant. In general, the parton shower cannot provide this, 
which means there will be even a 2 L^ ut ~ y / ftl terms induced in eq. (2.22). 

Pragmatically, the inclusive cross section can be restored to formal 0(a s ) accuracy by 
either explicitly including the corresponding a 2 corrections in the cumulant to cancel these 
terms, where numerical methods to do so have been described very recently in Refs. [29, 36], 



6 One can see this for example in the case of 2-jettiness in figure 3 in section 3. Here, the NLL'+L03 result 
is much closer to the slightly higher NNLL'+NL03 best prediction than the fixed NLO3 result. We have 
checked that in this case, adding the NLO3 to the NLL'+LC>3 by expanding it to 0(o? a ) forces the result to 
move in the wrong direction toward the lower NLO3. 
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or alternatively by explicitly restricting the fixed NLO corrections in the spectrum to the 
fixed-order region at large r, such that the induced O(a^) terms in the total cross section are 
not logarithmically enhanced and are formally 0(a 2 s ). This is essentially the approach taken 
in Ref [28]. However, neither of these approaches improves the perturbative accuracy in the 
spectrum outside the fixed-order region. 

2.2 What Monte Carlo Can Do for Resummation 

For T being the resolution variable between N and more than N jets, we showed in the 
previous subsection that combining the NLO matrix-element corrections for N and N + 1 
partons at the level of the singly differential T spectrum is equivalent to combining the NNLL 
resummation of the singular contributions with the higher-order nonsingular contributions. 
Our goal now is to extend this singly differential description to the fully differential case, in 
order to use the full iV-parton and (N + l)-parton information of the matrix elements. We 
will use the notation (N)LOjv or (N)LOn + i to indicate up to which fixed order in a s the 
./V-parton or (N + l)-parton matrix elements are included. 

To start with, it is straightforward to generalize the jet resolution spectrum da/dT and 
its cumulant o"(T" cut ) to include the full dependence on the iV-body Born phase space, 

da da 



dT d<S> N dT 
a(T-) - ^T^) = fdTj° df e(T<T^), (2.23) 



such that eq. (2.4) becomes 



d(7incl da (T cut ) + [dT d<T = 9(T > T cut ) . (2.24) 



d<I>Ar d$Ar J d® N dT 

Here, d<7i nc i/d$jv is the inclusive iV-jet cross section. The discussion in section 2.1 can be 
precisely repeated in this case, since the perturbative structure of the differential spectrum 
dcr/d$ArdT with respect to T is precisely the same as in eqs. (2.5) and (2.6). Namely, we can 
write it as the sum of singular and nonsingular contributions, 

da _ da sin g | da nons 
d$ N dT ~ d$ N dT + d^dT ' ( ' ' 

The nonsingular contributions are general functions of $tv and T, but as before are integrable 
in T for T — > 0. The singular contributions have the structure 



da sin s da B 



d® N dT d$jv 



CLi($jv, a s ) 6(T) + c n($N, a,) ^C n 0- 

n>0 y W 



(2.26) 



where das/d^M is now the fully differential Born cross section. Since the singular contribu- 
tions arise from the cancellation between virtual and real IR singularities, which only know 
about ^iv, their T dependence naturally factorizes from the $at kinematics of the underlying 
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hard process. This is what allows the resummation of the singular terms to higher orders 
for a given point in At LL, the entire <3?at dependence is that of the Born cross section. 
At higher logarithmic orders this is not the case anymore, since the coefficients C n can have 
nontrivial <3?at dependence. In addition, the precise definition of T also becomes important. 
Depending on its definition, the higher-order singular coefficients can depend on clustering 
effects or other types of nonglobal logarithms [37-42] which can be difficult to resum to high 
enough order with currently available methods. Therefore, it is important to choose a res- 
olution variable with simple resummation properties. An example is A-jettiness for which 
the complete NNLL resummation for arbitrary A is known [43, 44]. For the purpose of our 
discussion below, we will assume that a resummed result for the spectrum and its cumulant 
in eq. (2.24) at sufficiently high order is available to us. 

We can think of the cumulant dcr/d$Ar(7^ut) in eq. (2.24) as the exclusive A-jet cross sec- 
tion with no additional emissions (jets) above the scale T" cut , while the spectrum d<7/d<3?7vd7~ 
for T > 7~ cut is the corresponding inclusive (A + l)-jet cross section. While the cumulant 
da/d&N(T cut ) is differential in d$Ar, and thus already as differential as it can be, the spec- 
trum contains a projection from the full d<l>>Ar + i phase space down to d^ArdT. To also be 
fully differential in the (AT + l)-jet phase space, we can generalize eq. (2.24) to 



where da/d&N+i(T) denotes the fully differential spectrum for a given T = T(3>jv+i)- We 
explicitly denote the dependence on T and 7~ cut to clearly distinguish the spectrum from the 
cumulant. We have also used the short-hand notation 



where $at(<I>7v+i) denotes a projection from an (N + l)-body phase space point to an A- 
body phase space point. This projection defines what we mean by N jets at higher orders in 
perturbation theory. Note that beyond LO, both the cumulant da/d&]y(T cnt ) and spectrum 
da/d&N+i{T) must be well-defined jet cross sections, i.e., they require a specific IR-safe 
projection from to for both k = A and k = N + 1. We will see below where this 

definition enters. Using eq. (2.27) at the next-higher perturbative accuracy requires us to 
combine the higher-order resummation in T for the cumulant and spectrum with the fully 
exclusive A-jet and (A + l)-jet fixed-order calculations at NLOjv and NLOjv+i- To achieve 
this we have to construct appropriate expressions for the cumulant, dcr/d$Ar(T cut ), and the 
spectrum, dcr/d<I>jv+i(T), which we do in the next two subsections. 

2.2.1 Matched Cumulant 

We start by discussing the cumulant in eq. (2.27). Since the resummation is naturally differ- 
ential in the d$jv of the underlying Born process, we can combine the resummed result with 




(2.27) 



dj^v+i 
d$Af 



d$7V + i<5(3>Ar - <S> N (<S> N+1 )) 



(2.28) 
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the fixed-order one by adding the fixed-order nonsingular contributions to it, 



j j resum 

da (r cut ) = da t 
dq?Ar d<Pjv 



(2.29) 



FOJ 

The first term contains the resummed contributions, while the difference of the two terms in 
square brackets provides the remaining nonsingular corrections that have not already been 
included in the resummation. The NLOat fixed-order result is given by 

<lfT " NL ir* A ) = B N (* N ) + V N (* N ) + JdTe(T<T cut ) J^^ B N+i^N+i), (2.30) 



d$>Ar 

where Sjy, B^+i are the iV-parton and (N + l)-parton tree-level (Born) contributions, Vjv is 
the iV-parton one-loop virtual correction, and we abbreviated 

?f*tl = d^TV+i 6[T - T($n+i)} 6[$n - <S>n($n + i)} ■ (2.31) 
d<l?Ard/ 

Here, T(&n+i) implements the definition of T. The NLOat result also depends on the 
projection from <&n+i to <3?at, i.e., the precise NLO definition of $jv- However, this dependence 
only appears in the nonsingular corrections. For a given definition of T, the singular NLO 
corrections do not depend on how the remaining $>n+i phase space is projected onto $at, 
since they arise from the IR limit in which all (IR-safe) definitions agree. In eq. (2.29), the 
singular contributions inside the full fixed-order cumulant, d<7 FO /d < 3?jv(7 xu1: ) are canceled by 
the NLO expansion of the resummed result at NLL^. or higher, leaving only the nonsingular 
fixed-order contributions in square brackets. 

2.2.2 Matched Spectrum 

To properly combine the higher-order resummation in T with the fully differential (N + l)-jet 
fixed-order calculation, the inclusive (N + l)-jet spectrum dcr/d<I>Ar + i(T) in eq. (2.27) has to 
fulfill two basic matching conditions, 

Condition 1: ( — ^— (T) = , (2.32) 



Condition 2: — - (T) 



d$ 



7V+1 



_ da FO 

FO d<l>Ar + i 



(2.33) 



The first condition states that integrating the fully differential spectrum over the additional 
radiative phase space has to reproduce the correct spectrum in T including the desired resum- 
mation and fixed-order nonsingular corrections, such that eq. (2.27) reproduces eq. (2.24). 
The second condition states that the fixed-order expansion of the fully differential spectrum 
has to reproduce the full (N + l)-jet fixed-order calculation, where at NLOtv+i 

— - = B N+1 ($ N+l ) + V N+l ($ N+ i) + / -/i 2 B N+2 ($ N+2 ) . (2.34) 



- 16 - 



Here, -Bat+i and Bn+2 are the (N + l)-parton and (JV + 2)-parton tree-level (Born) contri- 
butions, and Vn+i is the (N + l)-parton one-loop virtual correction. Integrating over d^jy + 2 
in the last term now requires a projection from &n+2 to ^jv+i, 



d<£jv + 2 

d$7V+i 



d$ 



N+2 



6\$ 



N+l 



<S>N+l($N+2)] 



(2.35) 



analogous to eq. (2.28), which now defines precisely what we mean by N + 1 jets at NLO. 

In principle, there is some freedom to construct an expression for dcr/d<I>Ar+i(T) that 
satisfies both conditions to the order one is working. Our master formula to combine the 
resummed spectrum dcr resum /d$ArdT with the fully differential dcr FO /d$>N+\ is given by 



da 



d$jv+i 



(T) 



da 



FO 



d$> N +i 



d<7 r 



da r 



d$ N dT / d$ N dT 



FO. 



(2.36) 



Expanding the right-hand side to a given fixed order, we can see immediately that Condition 
2 is satisfied by construction. Imposing Condition 1 yields the consistency (or "matching") 
condition 



da 



d<f> N dT 



da 



FO 



da T 



da T 



d$ N dT/ d$ N dT 1 ' T ' 1T ' 

If the resummed result already has the nonsingular contributions at the desired fixed order 
added in, then the term in brackets is by construction equal to unity for any value of T '■ 
Otherwise, the expansion of the resummed result reproduces the singular terms of the full 
fixed-order result, leaving the nonsingular fixed-order contributions, such that we get 



da ( J cr sing,resum d<7 nons 

d$ N dT " d$ N dT + d<S> N dT 



do-sing.resum / (-[^sing 

d3>ArdT / d<S> N dT 



(2.38) 



Here d,j sm g> resum denotes the pure resummed result only containing the resummation of the 
singular contributions. Hence, eq. (2.36) not only multiplies in the additional dependence on 
&n+i/&n at fixed order, but if needed also adds the nonsingular corrections to the spectrum 
multiplied by the higher-order resummation factor. (Note that for the expansion of the 
resummed result to indeed reproduce all the singular terms at the desired fixed order, the 
resummation has to be carried out to sufficiently high order, which we have already seen in 
section 2.1.) 

To apply Condition 1 we have to integrate eq. (2.34) using the projection onto $tv and 
T in eq. (2.31). Therefore, to get the correct T spectrum at NLOjv+i, the projection in 
eq. (2.35) has to satisfy 

T[<S> N +1 {<S> N+2)] = T{§ N+ 2) , (2.39) 

i.e., it has to preserve the value of T when constructing the projected <&_/v+i point. Usually, the 
simplest way to handle this would be to use the left-hand side to define T(&n+2)- However, 
in our case eq. (2.39) provides a very nontrivial condition on the projection, since T(&n+2) 
is already defined by our choice of jet resolution variable, which in particular has to be 
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inclusive iV-jet 


exclusive iV-jet 


inclusive (N + l)-jet 


notation 


fixed, order 


accuracy 


log. order 


accuracy 


fixed order 


accuracy 


1jLjT-|-LjO 7Vr_l_1 




~ 1 


LL 


-1/2 

LXg 


LOa/"+i 


~ 1 


NLL r 


LOat 


~ 1 


NLL 


~ 1 






NLL r +LOAr + i 


LOat 


~ 1 


NLL 


~ 1 


LOa/+i 


~ 1 


NLL^+LOtv+i 


NLO^v 


~ a s 


NLL' 


1/2 


LOa/+i 


~ 1 


NNLL7-+NLO n+1 


NLO N 


~ a s 


NNLL 


~ a s 


NLOjv+i 


~ a-s 


NNLL^+NLOat+i 


NLO N 


~ a s 


NNLL' 


3/2 


NLOjv+i 





Table 1. Fixed and resummation orders and their achieved accuracy in a s . 



resummable. This turns out to be a nontrivial technical challenge one has to overcome to be 
able to satisfy Condition 1. We will see where this enters in section 3.1.2 and section 4.1.2. 

Note that to ensure that the resummation factor in square brackets in eq. (2.36) is well 
behaved in the fixed-order region at large T, it is important to turn off the resummation such 
that the ratio of the resummed spectrum and its expansion becomes 0(1) up to higher fixed- 
order correction. In principle, the fixed-order result in the denominator can also become 
negative at very small values of T . This is not a problem in practice, since this region is 
explicitly avoided by imposing the cut T > T cnt - 

2.2.3 Perturbative Accuracy and Order Counting 

The appropriate order counting in the resummation and fixed-order regions is precisely the 
same as in section 2.1, so there is no need to repeat it here. Applying eq. (2.36) at the very 
lowest order, namely LL CT resummation with LOat+i fixed-order corrections, we get 



da 



>N+1 



d$ 



N+l 



r>o 



B N+1 ($ N+1 ) exp[^a s L 2 ] , (2.40) 



where Bn + i($>n+i) scales as a s /T relative to Bn(&n) at small T, and we used that at LL CT 
the ratio in brackets in eq. (2.36) is just the resummation exponent. This directly corresponds 
to the CKKW-L procedure [1, 2, 33, 34], which multiplies the tree-level matrix elements with 
the shower Sudakov factors. Hence, we can think of our master formula eq. (2.36) as a 
consistent extension of this to higher orders. 

As demonstrated in section 2.1, going to the next higher perturbative accuracy in all 
phase space regions requires the NLL^. and NNLL CT series of logarithms. We obtain these by 
performing the full NLL' and NNLL resummation in the exponent, as well as the fixed NLOat 
and NLOat+i corrections in the cumulant and spectrum, respectively. The resummation 
naturally connects both jet multiplicities, since the NLOat corrections are included in the 
cumulant and are part of the resummation for the spectrum starting at NLL', where they 
effectively predict the singular NLOat+i contributions, and the full NLOat+i corrections are 
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obtained by adding the nonsingular corrections to the spectrum. In the following, we will 
use the notation (N)NLL^+(N)LOat+i to indicate the resummation order for the employed 
jet resolution variable together with the (N + l)-jet fixed order. For simplicity we do not 
explicitly denote the A-jet fixed order and keep it implicit in the resummation order, i.e., 
LOtv at (N)LL and NLOjv at NLL' and above. This is summarized in Table 1. 

An immediate and important question to ask is to what accuracy resummed spectra for 
jet resolution variables other than T are predicted in our approach. A detailed theoretical 
investigation of the formal resummation order one attains for other variables would be very 
interesting, but is beyond the scope of the present work. What is certainly clear is that 
other variables will not be resummed at the same formal level as the primary jet resolution 
variable T itself. However, we know that other variables are correct to NLOat+i, while at the 
same time the inclusive cross section is not changed, as it is independent of which variable 
one integrates over. This implies that the NLOyv+i corrections for other variables do not 
induce uncanceled higher-order logarithmic terms as in eq. (2.22), and hence some higher- 
order resummation must be partially retained for other observables as well. Numerically, the 
higher-order resummation in T provides an improved weighting of the IR region of phase 
space, from which other variables are expected to benefit as well. We can validate to what 
accuracy other variables are obtained by comparing predictions from our highest order to the 
analytically resummed results for other observables, which we do in section 3.3. 

2.3 Extension to More Jet Multiplicities 

The method proposed in this paper is completely general and can be extended to more jet 
multiplicities, essentially by iterating the procedure discussed in section 2.2. We start by 
introducing separate jet resolution variables 7jv to distinguish A vs A + 1 jets, 7jv+i to 
distinguish A + 1 vs A + 2 jets, and so on. One can choose any IR-safe observable that 
goes to zero in the limit of A pencil-like jets. For each A the inclusive A-jet cross section is 
obtained by combining the cumulant and spectrum for 7~n as in eq. (2.27), 

dcTind dcr ,_. cutx , f d$ N+1 da 



J d<Pjv d$Ar + i 



d$ N d<3?jv 
do "mcl da t [ d® N+2 da 

Hcb = Hcb Vn+1) + / ^ ^ {In+i)V{In+1 > l N+ D , 

d®N+l d<I>Ar + i J d<PAT + i d3>7v+2 



dOincl dcr ^ ( . u| 



-> oo) • (2.41) 



The exception is the highest jet multiplicity, A max , for which T^ ut = oo, corresponding to 
the fact that no additional jets are resolved. 

For the cumulants in eq. (2.41), the discussion in section 2.2.1 applies separately for each 
N, so the cumulants matched to higher resummed and fixed order are given as in eq. (2.29) 
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by 



da 



N 



(r N ut ) 



da r 



d$ 



N 



+ 



da 



FO 



d$ 



A? 



-(rr) 



da r 



d$ 



v 



-(7? 



•cut 1 
N / 



FO. 



(2.42) 



The fully differential T/v spectra d<7/d$jv+i(7w) are now obtained recursively as follows. 
We start with the highest jet multiplicity, iV maX) for which no resummation is needed since 
7jV ^ s essentially removed. Furthermore, the highest jet multiplicity is by construction only 
required at leading order, where the result is simply given by the Born contribution, 



da 



-(7J 



cut 



OO 



d$ 



N n 



Bat ($tv ) 

1 v max V 1 v max / 



(2.43) 



For each N < iV max we apply the discussion in section 2.2.2. To combine the resumma- 
tion in T/v with the (N + l)-jet fixed-order calculation, the fully differential T/v spectrum 
d<r/d<l>/v-|_i(Tv) must satisfy the matching conditions as in eqs. (2.32) and (2.33), 
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These can be satisfied by a straightforward generalization of eq. (2.36), 
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(2.44) 
(2.45) 

(2.46) 



The prefactor on the right-hand side is now the inclusive (N + l)-jet cross section from 
eq. (2.41). This is what ties together the different jet multiplicities. The condition in eq. (2.45) 
now leads to the consistency condition 
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(2.47) 



which states that for each N the cumulant and spectrum in T/v+i must be included to suffi- 
ciently high order such as to reproduce the (N + l)-jet fixed order that is required by the T/v 
spectrum. Imposing the condition in eq. (2.44) yields the consistency condition for the T/v 
spectrum 
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(2.48) 



d$/vdTv d$>/v+i/ d$/vdT/v 
which is the generalization of eq. (2.37). To satisfy eq. (2.47) at NLO/v+i it requires that 

T N {$N+l($N+2)} = Tv($/v+ 2 ) , (2.49) 

as in eq. (2.39). That is, for each N the projection from <&/v+2 to <&/v+i which defines the 
(N + l)-jet cross section at NLO has to preserve the value of T/v- In addition, eq. (2.48) 
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requires that upon integration the Tn+i resummation contained in do"i nc i/d<I>jv + i does not 
interfere with the Tn resummation, e.g., by inducing higher-order logarithms in T/v- Since 
eq. (2.41) relates the da mcl /d&N to dcr/d^jv+j^T/v), the relationship in eq. (2.46) gives rise 
to a recursive definition, which combined with the result for the highest jet multiplicity in 
eq. (2.43) determines d<r/d<I>Ar + i(7Ar) for all N. 

In the Monte Carlo implementation, the phase space is split up recursively as 

-Un )+ — -r^ 0( In > In ), 



d$N d$N J d$Af d$A? + i 

^cb = H* y'N+i) + / ^ V\Jn+\ > 7jv+i) , (2-50) 

d^jy+i a<I>Ar + i J d®N+i a<S N+2 



where in each step, the total cross section for N or more jets is separated into an exclusive 
iV-jet cross section, which is assigned to partonic events with N final-state partons, and the 
integral over the remaining cross section for N + l and more jets. For the highest multiplicity, 
iVmax, the remaining cross section for N max or more jets is represented by events with iV max 
final-state partons. 

Note that the structure of eq. (2.50) is very similar to eq. (2.41). The crucial difference 
is that in eq. (2.50) each inclusive cross section on the left-hand side is the same that appears 
under the integral on the right-hand side in the line above. By comparing eq. (2.50) with 
eq. (2.41) and repeatedly inserting eq. (2.46) we obtain the higher-order "fully-resummed" 
exclusive iV-jet cross sections that serve as inputs to the Monte Carlo. Abbreviating the 
resummation factor in eq. (2.46) as 



^^.rcsum / ^^.resum 
UN ^ N ' TN) = d^dfk/ d^dfk 



(2.51) 



FO 



we obtain 



rl,T MC H/t 

W " } ~d^ {TN h 

J MC j 



(7^1 oo) U n ($ n ,Tn) U N+1 ($ N+lt TN+ij 



x • • • x ^ m a X -l(^JV m a X -l,T7V max -l) . (2.52) 

The careful reader will have noticed that the above is in one-to-one correspondence to the 
structure generated by a parton shower with up to iV max emissions. The crucial difference is 
that in our case all ingredients are well-defined physical jet cross sections defined in terms 
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of a global jet resolution variable. This allows us to systematically increase the perturbative 
accuracy by computing the relevant ingredients to higher order in resummed and fixed-order 
perturbation theory as well as to systematically estimate the perturbative uncertainties. The 
analogous parton-shower-like structure underlies the CKKW-L ME/PS merging, which re- 
places the splitting functions in the shower with the full tree- level matrix elements. Restrict- 
ing eq. (2.52) to the lowest order as in eq. (2.40), it reduces to the ME/PS merging as a 
special case. 

In principle the above construction allows us go to even higher fixed and resummation 
order, as long as the fixed-order ingredients are available and the resummation is known to 
a correspondingly high enough order. It also lets us combine as many jet multiplicities as 
we like at the order they are available. In particular, it is straightforward to add additional 
multiplicities at the lowest accuracy in a CKKW-L-like fashion. 

2.4 Attaching Parton Showering and Hadronization 

In the Monte Carlo, a point in <3?jv is represented by N (massless) four- vectors together with 
the appropriate flavor information. We then generate events with N to -/V max partons and 
assign the iV-parton events the weight dcr^ c /d$Ar(T^ ut ), the (_/V + l)-parton events the weight 
dcr^ 1 /d < I ) Ar + i(T^ l ^ 1 ), and so on. The events with iV max partons are assigned the weight 
^>N m J^N n ^ = B Nmax (* Na J. The 6(T N > 7^ ut ) functions in eq. (2.50) are included 
in the weight, which means that all events with > N + 1 partons that have TJv < 7Jv ut get 
zero weight. 7 In this way, summing up the weights of all events we can integrate up the 
cross sections in eq. (2.50), including arbitrary kinematic cuts in <E>7v> &N+1, etc. What is 
important is that although the events contain massless partons, they represent the exclusive 
jet cross sections of eq. (2.52). (From the resummation point of view, the massless partons 
represent the kinematics of the hard function.) 

In the next step, the events are given as a starting point to a parton shower, whose 
purpose it is to fill up the jets with additional emissions inside the jets. Formally, this 
means that the shower should not be allowed to change the underlying distribution in the jet 
resolution variable, since this has already been computed at the higher perturbative accuracy. 
For example, starting from an event with N + 1 partons with kinematics # n+i and weight 
dc^ 1 /d ( I ) Ar + i(T^ u ^ 1 ), the fully showered event should have the same jet kinematics &n+i as 
the unshowered event from which it originated. Most importantly, the showered event should 
have the same value of Tn($n+i) and should have T/v+i < 7tv+i> so ^ s ^ m has the correct 
weight d<r)J-'^ 1 /d<I ) Ar_|_i(T^ U | t 1 ). (The shower itself does not change the weight of the event.) In 
practice, this is quite a nontrivial constraint on the shower. The easiest way to enforce this is 
to repeatedly run the shower on the same event until it produces an acceptable showered event, 
where we allow the value of T/v to be changed at most by a small numerical amount consistent 
with a power correction. This method is of course computationally intensive (though it is 

technically, the split up of phase space is usually flavor-aware. This means that an event with Tzv < T/v ut 
is only set to zero if the closest two partons produce a QCD singularity. 
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not computationally prohibitive), since one may have to rerun the shower many times, and 
it would be interesting to develop a more efficient way of constraining the shower for this 
purpose. Notice that in this procedure no events are discarded so the cross section is not 
changed. 

In the final step, the showered event is passed to the hadronization routine. In this case, 
there are no constraints on the kinematics of the hadronized event, i.e., the hadronization is 
allowed to smear out the T/v spectrum. The reason is that our perturbative calculation does 
not take into account nonperturbative effects, but are instead supplied by the hadronization. 
This is discussed in more detail in section 3.1.3. 

3 Application to e + e Collisions 

In this section, we apply the framework described in section 2 to e + e~ — ^ 2/3 jets, imple- 
mented in the Geneva Monte Carlo. The higher-order resummation for 2-jet event shapes in 
e + e~ collisions is very well understood and many precise measurements from LEP exist, which 
are used, for example, for precise determinations of the strong coupling constant a s [45-51]. 

In this context, one important aspect is the interplay between both resummed and fixed- 
order perturbative contributions with the nonperturbative corrections. Here, the Geneva 
framework provides an important development by being able to combine the perturbative 
higher-order resummation with the nonperturbative information provided by Pythia's hadroniza- 
tion model [52, 53]. For example, this allows us to use a common theoretical framework to 
make predictions for different phase space regions and different observables. 

The e + e~ implementation also provides an important and powerful validation of our 
approach and its practical feasibility, while avoiding the additional complications arising for 
hadronic collisions, such as initial-state radiation and parton distribution functions (PDFs). 
The implementation and first results for pp collisions are presented in section 4. 

In our e + e~ implementation, we use 2-jettiness, T2, as the 2-jet resolution variable, which 
is defined as [43] 



and is simply related to thrust T [30] by 7i = E cra (l — T). Its kinematic limits are < 7i < 
-E'cm/2. In the limit T2 — > 0, there are precisely 2 pencil-like jets in the final state, while for 
75 ~ E cm , there are 3 or more jets. We perform the resummation in 7i to NNLL' and include 
the full NLO2, NLO3, and LO4 fixed-order matrix elements, i.e., we obtain NNLL^+NL03 
predictions. 

The default running parameters for our e + e~ studies are E cra = 91.2 GeV, a s {mz) = 
0.1135, and Pythia 8.170 with e + e~ tune 1. Using this value of a s (mz) is motivated by 
the fact that it was obtained from fits to the thrust spectrum using N 3 LL/ resummation. 

8 The a s value used inside Pythia's parton shower is not changed from the value set in the tune. This 
is not inconsistent, since here the strong coupling functions as a phenomenological parameter, regulating the 
amount of showering. 





- 23 - 



These fits were performed in a region (corresponding to 6 GeV < 7i < 30 GeV for our .E cm ) 
where the nonperturbative corrections due to hadronization are power suppressed and can 
be described by a single nonperturbative parameter, which leads to a shift in the spectrum 
and is included in the fit in Ref. [49]. We find that this value of a s (mz), in conjunction with 
Pythia's tune 1, provides overall the best description of the data, including the peak region 
below 7i < 6 GeV and other 2-jet event shapes. For comparison, we show results using the 
world average a s (mz) = 0.1184 [54] as well as from using Pythia tune 3. 

In the next subsection, we summarize the various ingredients that go into the master 
formula, with the intention of giving a concise and informative overview, while leaving a de- 
tailed discussion of our implementation to a separate publication. In section 3.2, we discuss 
the 7i spectrum, validating our implementation using analytic predictions as well as compar- 
ing our results to LEP data. In section 3.3, we present our results for other 2-jet variables, 
namely C-parameter, heavy jet mass, and jet broadening, comparing Geneva's predictions 
at NNLL^+NL03 to the analytic higher-order resummation for each variable as well as to 
the experimental measurements. In all cases, we find good consistency and agreement with 
the data. 

3.1 Ingredients 

The master formula is given by 

d<7 inc i da t f d<J> 3 da ln -, aln - ^ _ cut 



where 



da ,„ — tN da r 



-C7? Ut ) = -^(7? ut ) 



d$2 d$2 

da daM A1t 1,s,u1 /da'' 

-{12) 



d$3 V ' d$ 3 \ d$2dT2/ d$ 2 dr 2 



(3.3) 

FO/ 



Its three key ingredients are the higher-order resummation of 2-jettiness, which we include 
at NNLU7— I-LO3, the full fixed-order matrix elements at NLO2, NLO3, and LO4, and the 
interface to parton showering and hadronization, for which we use Pythia 8. 

Following the construction in section 2.3 with iV max = 4, the inclusive 3-jet cross section 
is separated into 3 and 4 or more jet contributions using 3-jettiness, 7*j, as our 3-jet resolution 
variable, 

For e + e~ collisions, iV-jettiness is defined by [43] 

T N = min (-£fc - fn ■ Pk) , (3.5) 
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Fixed-order corrections Resummation input 





singular 


nonsingular 




Fcusp 




LL 


L0 2 






1-loop 


1-loop 


NLL 


L0 2 




1-loop 


2-loop 


2-loop 


NLL' 


NL0 2 




1-loop 


2-loop 


2-loop 


NLI/+LO3 


NL0 2 


L0 3 


1-loop 


2-loop 


2-loop 


NNLL+LO3 


NL0 2 


L0 3 


2-loop 


3-loop 


3-loop 


NNLL' 


NNLO2 




2-loop 


3-loop 


3-loop 


NNLI/+NLO3 


NNLO2 


NLO3 


2-loop 


3-loop 


3-loop 



Table 2. Perturbative inputs included at a given order in resummed and fixed-order perturbation 
theory. The columns in the resummation input refer to the noncusp anomalous dimension (7^), the 
cusp anomalous dimension (r cusp ), and the QCD beta function (/?). 



where i = 1, • • • ,N and hi is a unit vector along the direction of the ith jet, where the jet 
directions can be determined by a jet algorithm or by directly minimizing Tn- 9 There are N 
pencil-like jets in the limit Tn — > and N or more jets in the limit Tn ~ -^cm- 

As discussed in section 2.3, the master formula naturally incorporates the resummation of 
the 3-jet resolution variable in eq. (3.4) and extends to higher jet multiplicities, i.e., iV max > 4. 
However, since our current focus is on the main conceptual development of combining the 
higher-order resummation with the fixed NLO matrix elements for 2 and 3 jets, we leave these 
extensions to future work. As we will not be interested in the T3 spectrum or other exclusive 
3-jet observables, it is sufficient for our purposes to calculate the two terms on the right-hand 
side of eq. (3.4) at fixed order. Thus, we use 

^W^^V'). (3.6) 

In the results that follow, we use 73 CUt = 2 GeV and a T£ nt value between 0.5 — 1 GeV, which is 
selected randomly from a flat distribution. This smoothing out of 72 CUt avoids small numerical 
discontinuities that can arise with a sharp cutoff. 

3.1.1 Resummation 

Our jet resolution variable, 7i, has the important property that it can be factorized. The 
factorization theorem for the T2 spectrum provides the resummed prediction that is one of the 
primary inputs to our master formula in eq. (3.3). It is obtained by using the framework of Soft 
Collinear Effective Theory (SCET) [57-60] and allows the resummation to be systematically 



This definition agrees with T2 in eq. (3.1) for massless final state particles, which is the limit in which 
resummation is carried out. It does affect the nonperturbative corrections when including hadron masses [55, 
56]. We use the definition of T2 in eq. (3.1) to be able to directly compare to the experimental data for thrust. 



- 25 - 



carried out to higher orders and combined with the nonsingular fixed-order result. Our 
highest-order resummed input to the master formula has NNLL' resummation. We use the 
standard resummation formalism, where the large logarithms are resummed in the exponent 
of the cross section, with the corresponding resummation orders summarized in table 2. 
We write the jet resolution distribution in 7i as 

dcrf sum _ dcr^ ing t dcrf ns 

d$ 2 dT 2 ~ dn 2 dT 2 + dQ 2 dT 2 ' ( ' J 

where the separation into singular and nonsingular contributions was discussed in section 2.2 
(see eq. (2.25)). The singular contribution is given by [61, 62] 

d0"2 1Ilg diJB /rn O f T I \ T / n (rr s i s 2 



dn 2 dT 2 dQ 2 



H 2 (Ec m ,n) dsids 2 Ji(s 1 ,fi) J 2 (s 2 ,fi) S 2 (T 2 - (3.8) 



Here d<& 2 = d£l 2 = d cos 6d(j) is the angular phase space for the orientation of the thrust axis 
with respect to the beam and doe/df^ is the tree- level 2-parton cross section. Note that 
the overall dependence on £l 2 here is that of the Born cross section, which is correct in the 
limit T 2 — > in which the eq. (3.8) is obtained. The hard function H 2 in eq. (3.8) contains 
the fixed-order 2-parton matrix elements, which describe the short-distance corrections at 
the scale E cm . The jet functions J\ and J 2 describe the back-to-back collinear final-state 
radiation along the thrust axis and the soft function S 2 describes the soft radiation between 
the jets. The soft function contains perturbative and nonperturbative components, which can 
be separated as [63-65] 

S 2 (T 2 ,fi) = JdkSr\T 2 -k^)f{k^) } (3.9) 

where Sf 6 (7a — k,fx) is the perturbative soft function, while the shape function f(k,fi) de- 
scribes the nonperturbative hadronization corrections. For T 2 ~ Aqcdj the shape function 
gives an 0(1) contribution to the cross section, while for T 2 3> Aqcd> it can be expanded 
and only the leading 0(Aqcd/72) nonperturbative power correction is relevant. For further 
discussion and the derivation of the factorization theorem, see Refs. [61, 62]. The resummed 
prediction used in Geneva only includes the perturbative soft function, while the nonpertur- 
bative corrections are provided by the hadronization in Pythia. 

The nonsingular contribution in eq. (3.7) is given by the spectrum at fixed order with the 
singular terms subtracted. It includes all 0(T 2 / E cra ) corrections to the singular distribution 
to a given order in a s . The 0(a s ) nonsingular corrections in T 2 are known analytically and can 
be taken from Ref. [49], so we include them in our resummed result. Each function in eq. (3.8) 
depends on the renormalization scale [i and the characteristic scale of the physics it describes. 
These are [iu ~ -Ecm, fJ-j ~ \/T 2 E cmi and fis ~ T 2 for the hard, jet, and soft functions, 
respectively. Renormalization group evolution (RGE) between the soft, collinear, and hard 
scales resums the logarithms of the form Ihhs/^h ~ ^T 2 /E cm and ln/ij///^ ~ lnT 2 /E cra 



-26- 



in the factorized singular distribution in eq. (3.8). The anomalous dimensions and singular 
fixed-order corrections required at a given resummation order are summarized in table 2. 

The resummed cumulant in eq. (3.3) is obtained in an analogous way to the resummed 
T~2 distribution. It is given by a singular and nonsingular component 



where the singular contribution is obtained by integrating eq. (3.8) over 7i from to 72 CUt - The 
nonsingular contribution to the cumulant is given by the difference between the fixed-order 
result and the resummed singular terms expanded to fixed order. 

The perturbative uncertainties in the resummed spectrum are estimated by scale vari- 
ation and receive a contribution from two distinct sources, the fixed-order corrections and 
the higher-order logarithmic resummation. The fixed-order uncertainties are estimated by a 
correlated overall variation of all scales by factors of two. The resummation uncertainties are 
instead estimated by varying the lower scales \xj^Ti) and ns^Vi)-, which are functions of 72, 
and are referred to as profile scales [49, 65, 66]. The profile scales satisfy the criteria that 
in the resummation region, fJ^jsiTi) have their canonical scaling (given above) and in the 
fixed-order region, fJ>j t s{T2) ~ HH, which turns off the resummation. In the transition region 
the profile scales provide a smooth interpolation between the resummation and fixed-order 
regions. These three regions are determined based on where the fixed-order singular contribu- 
tions dominate over the nonsingular ones. The variations in the profile scales subject to the 
above constraints determine the resummation uncertainty, where we take the largest absolute 
variation from the central scale. The resummation uncertainties are combined in quadrature 
with the fixed-order uncertainties to generate our theory uncertainty estimate. For a given 
partonic event in Geneva each profile scale variation gives rise to a different event weight, 
which is computed analytically. Hence, we can provide each event with its own perturbative 
uncertainty estimate by assigning it several weights from the profile scale variation in addition 
to its central weight. 

3.1.2 Fixed Order 

As we can see from eqs. (3.4) and (3.6), we need the 3-jet cumulant da/d^s(T^ ut ), as well 
as the Born 4-parton cross section i?4(<&4). The Born 4-parton cross section is trivial and 
requires no further discussion. To calculate the 3-jet cumulant at NLO3, we use the generic 
formula given in eq. (2.34) 




(3.10) 





where 



d$4 



d$4^[*3 - $3 



(*4)] • 




d$ 3 
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The projection 3>g ($4) defines what we mean by $3 at NLO3. It implicitly depends on our 
choice of resolution variable, since eq. (2.39) requires it to satisfy 

7^[* 3 r ($ 4 )]=75(*4). (3.13) 

To deal with the IR singularities that are present both in V3 and in the integral of B4 
over $4, we use the FKS subtraction method [67]. We introduce a set of projecting functions, 
#m( < l ) 4), that partition the phase space into nonoverlapping regions, such that £ m 0m( $ 4) = 1. 
In our case this partition is effectively determined by the resolution variable, which is indicated 
by the superscript. The resulting partition must be such that each region m contains at most 
one collinear and one soft singularity. Then, we can write 

d<D 4 = Y, d$ 3 d$™ d 0£($ 3 , <*>rad) , (3.14) 

m 

where $™ d denotes the radiative phase space describing a 1 — > 2 splitting in each region. 

For each region m we define a mapping that identifies which particle in $3 is undergoing 
the 1—^2 splitting which generates the 

$Z = <S>T($3,^ d ) (3-15) 

phase space point. It also unequivocally defines how the recoil is shared amongst the remaining 
particles in the event, which is needed to enforce total momentum conservation. Notice that 
our definition of the ^-functions leaves us the freedom to include phase space regions in the 
partition which do not contain any IR singularity. This freedom is in fact essential to be able 
to satisfy Condition 1 in eq. (2.32), namely, to ensure that in each phase space region the 
correct functional form for the resolution variable 72 ($4) is used. For example, in e + e~ —> qqg 
production, the region in which the quark and antiquark are closest to each other does not 
contain any QCD singularity. Nevertheless, it must be treated as a separate region in the 
phase space partition, since in this region the invariant mass between the quark and antiquark 
determines the value of 7i- 

With this notation, one can write 



da 
d$3 



B 3 (<t> 3 ) + %(* S ) + ]T / d<S>™ d B 4 ($T)0l(<S> 3 , <O)0(T 3 < 7? ut ) • (3-16) 



If the region m contains an IR divergence, the FKS subtraction requires one to define the 
soft, collinear, and soft-collinear limits of which we denote as $^ d , ^^d> and ^racT' 

respectively, together with the resulting points in the 4-body phase space $™' s , &™' c , and 

n 
4 

da 



&T' CS - We can then write 



(3.17) 



d$ ms d$ mc d$ m \ CS 

d*m^(*4 )^m( 3) ^ ra d) d$^j > + ^ 
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where 6^ encodes the soft limit of the (^-functions and we have used the fact that in the 
collinear and soft-collinear limits the ^-functions are trivially satisfied. Also, since in each 
of these limits T3 = 0, the #(75 < 7^ cut ) functions are satisfied by construction. 

If m is not singular, in principle no such subtraction is needed and one could simply 
evaluate the 4-parton tree-level matrix element Z?4($™). However, given that the integral 
of the subtraction counterterms over the whole phase space is known analytically for both 
massless and massive partons [17, 67, 68] 

WO =E [/ Cd^WICl^,^) + / dCci c ^($D 

'd^f£ 4 ($r s )l> ( 3 - 18 ) 



we found it easier not to restrict the integration of the subtraction counterterms only in the 
singular regions of phase space but to extend it across all of phase space. 10 This ensures the 
complete cancellation of and the independence of the final results from the subtraction terms. 
The procedure outlined above takes care of all IR divergences, making the integrand in the 
square brackets of eq. (3.17), as well as the sum of V(<&s) + 1(^3), IR finite. 

The crucial point, discussed in section 2.2.2, is that our construction requires the phase 
space map that generates $>™($3, to preserve the value of 75, i.e. 

r 2 [^(*3,*Sd)] = r 2 (*3). (3.19) 

Comparing this to eq. (3.13) we see that the map $>™($>3, $™ d ) must be precisely the inverse 
of <I»3 _ (<I>4) in the region m. In principle, this condition can be relaxed to only hold up to 
power corrections. Additionally, the map can fail to preserve 75 in a region of phase space 
that gives a power suppressed contribution to the cross section. The phase space maps used 
in the standard FKS implementations [16, 67] were not designed to preserve the value of 
75, and thus they change its value by an 0(75 /75) amount over a large region of phase 
space 11 . Since 4-parton events with 75 < 73 CUt and 75 > are the only real emission 

contributions included in the NLO calculation for the 3-jet cumulant, eq. (3.11), one can 
impose the restriction 75 cut <C TJ" 1 and use the standard FKS phase space maps. However, 
this hierarchy strongly restricts 75 cut , and it is preferable to define a map that is specifically 
designed for our goals. We have constructed such a map, which preserves the exact value 
of 75 up to power corrections, except in a region of phase space whose contribution to the 
cross section scales as 0(T2/E cm ). In this region, the value of 75 is altered by an O (75/75) 
amount, meaning the net correction scales as 0{Tz/ E cm ). Therefore, enforcing the much 
looser constraint T^ nt <C E cm is sufficient to achieve our purposes. We postpone the detailed 
discussion of this map to a dedicated publication describing the implementation of Geneva. 



10 These integrals can also easily be denned by restricting the integration on the FKS variable £ up to some 
£cut value. These are however not in direct correspondence to the partition of phase space we are considering. 

11 Generically, an emission that takes a 3-parton event to a 4-parton event will change T2 by the scaling 
T2(4>4)-T 2 (4>3) ~T 3 . 
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3.1.3 Parton Shower and Hadronization 

The phase space points $2 and $3 represent jet kinematics, which are defined by the jet 
resolution variable, 2-jettiness. As discussed in section 2.4, we require that the parton shower 
does not change the underlying hard jet distribution so that the higher-order weights we 
calculate, da / d^2(T 2 cut ) and dcr/d$3(72), are correctly assigned. The parton shower, however, 
does not preserve the value of 7i- We address this problem in our current implementation 
in a physically motivated way. For small 7i the resummed singular jet resolution spectrum 
dominates and is determined up to power corrections of order A ~ T2/E cm . We require that 
for $>3 events, the change in T2 due to showering, AT2, satisfies AT2/T2 < A. This represents 
a power correction to the 2-jettiness spectrum, which scales as 1/72 for small 72- We also 
require that $2 events, which have 72 = when unshowered, remain in the 2-jet bin after 
showering up to a power suppressed correction, with 7i < 72 CUt (l + A'). The results shown in 
this section use A = 2 A' = 0.05. 

Furthermore, the shower must also be restricted to not change the NLO3 result. This 
requires that for 3-parton events we effectively only allow showers to start from 73° ut . Likewise 
we limit the showering of 4-partons events down from their 73(^4) value. This can be seen 
as a proxy for what would be the correct approach in a 7/v-ordered shower. 

We use Pythia 8.170 with e + e~ tune 1 for showering and hadronization. The choice of 
tune for e + e~ data in Pythia affects both the time-like showering and hadronization model. 
However, since in our implementation we restrict the shower from changing the T2 spectrum, 
the effect of changing the tune in Pythia primarily reflects the uncertainty from hadroniza- 
tion in Geneva. We have checked that this is also the case for observables other than 7i, by 
verifying that the effect of the tune on the showered Geneva predictions is very small com- 
pared to the change due to hadronization. The uncertainty from hadronization is associated 
with the nonperturbative contribution to the soft function in eq. (3.9) in our framework and 
is not included in our event-by-event perturbative uncertainties. As an indication of the size 
of the uncertainty from hadronization we also show Geneva hadronized results using e + e~ 
tune 3. 

It is important to note that we use the standard tunes in Pythia, without changing any 
internal parameters. Since in our approach the shower evolution in standalone Pythia 8 is 
substituted with higher-order resummation above the 2-jet resolution scale, we advocate that 
a separate tuning of Geneva + Pythia 8 should be employed to obtain the best results. 

3.2 Validation Using the Jet Resolution Spectrum 

Before comparing the Geneva prediction for various e + e~ spectra to analytic predictions 
and LEP data, we first validate the implementation of our procedure to combine higher-order 
resummation and full NLO matrix-element corrections by using the jet resolution spectrum. 
At the level of the singly differential 7i spectrum only, the standard approach to resummation 
achieves the same matching between resummation and fixed order by adding the nonsingular 
contribution to the resummed result. This provides a nontrivial cross check of the mas- 



- 30 - 



ter formula and in particular validates the event-by-event theory uncertainties generated by 
Geneva. For each comparison in this section, we show the peak, transition, and, tail regions, 
described in section 2.1, at the LEP center-of-mass energy E cm = 91.2 GeV. In all cases the 
error bars or bands on the Geneva histograms are built from its event-by-event perturbative 
uncertainties. The statistical uncertainties from Monte-Carlo integration are much smaller 
and not shown. 

3.2.1 Partonic Results 

The analytic resummed 72 spectrum is shown in figure 2 at successively higher orders: NLL, 
NLL/+LO3, and NNLL'+NL03 (see table 2 for the order counting definitions). The per- 
turbative uncertainties are generated by using the same profile scale variations employed in 
Geneva and discussed in section 3.1.1. The theory uncertainties decrease at increasing order 
and demonstrate excellent convergence at all values of 75- Below 7i < 0.5 GeV, we enter 
the purely nonperturbative region and the scale uncertainties diverge since even resummed 
perturbation theory breaks down. In the far tail the scale uncertainties also grow rapidly, 
which reflects missing higher fixed-order corrections. The uncertainties in the NNLL/+NLO3 
prediction diverge past the 3 parton endpoint at E cm /3, where the fixed-order prediction 
is only correct at leading order for 4 partons. In the transition region there is a smooth 
interpolation between the resummation and fixed-order regions. 

In figure 3 we compare the partonic 7i spectrum from Geneva with 7~2 Ut = 1 GeV to the 
analytic resummed results from figure 2. To illustrate the interpolation between resummed 
and fixed-order results we also show the pure resummed results at NLL' and NNLL' and the 
pure fixed-order contribution at LO3 and NLO3. The latter are calculated using Event2 [69, 
70] which serves as an independent cross check of our NLO3 implementation. Using the 
NLL'+L03 resummation of Ti and the LO3 fixed-order contribution as inputs to our master 
formula for the spectrum in eq. (3.3), the d(T FO /d$ 3 and da resum /d$ 2 d7i| FO contributions 
exactly cancel for the 7i spectrum. As a result, we see precise agreement between Geneva 
and the analytic NLL'+L03 result in figure 3(a)-3(c) in the peak, transition, and tail regions. 
This result agrees well with the pure NLL' resummed contribution in the peak, while in 
the tail it is consistent within uncertainties with the LO3 result, where the latter clearly 
underestimates the full perturbative uncertainties. 

At next higher order using as inputs to the master formula the NNLL'+L03 resummation 
of 7i and the NLO3 fixed-order calculation, we see that the central value and event-by-event 
uncertainties in Geneva agree very well with the full analytic NNLL'+NL03 resummed 
prediction in the peak and transition regions, as shown in figures 3(d) and 3(e). In the tail 
region, figure 3(f), Geneva has significantly smaller uncertainties of the same size as the pure 
fixed-order contribution. This is because there is a substantial cancellation between singular 
and nonsingular contributions in this region, which is incorporated differently in the analytic 
resummation and the master formula at NNLL'+NL03. For the former, the nonsingular a% 
contribution are added. This preserves the absolute size of residual resummation uncertain- 
ties, which are very small relative to the singular contributions but large relative to the total 
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Figure 2. Analytic resummation of T2 matched to fixed order. The central value is shown along 
with the band from scale uncertainties, as discussed in section 3.1.1, at NLL, NLL'+LC>3, and 
NNLI/+NLO3. 



result after cancellation. In the master formula in eq. (3.3), the nonsingular contributions are 
incorporated multiplicatively through the ratio of dcr FO /d$3 and d<7 resum /d3>2d72| FO - This 
preserves the relative size of residual resummation uncertainties, thus leading to much smaller 
absolute variations when compared to the final result. Comparing the Geneva prediction 
with the pure NNLL' resummed and NLO3 fixed-order results, we see that the master formula 
precisely interpolates as expected between the fixed-order and resummation regions, with the 
transition region properly describing the transition between the two, including uncertainties. 

Combining the exclusive 2-jet cross section with the integral of the inclusive 3-jet cross 
section, the Geneva prediction at NNLL/+NLO3 formally reproduces the total inclusive cross 
section at NLO. Numerically, we have cr^ot = 44.1 ± 0.2 nb. With T<z nt smeared between 
0.5 — 1 GeV the total inclusive cross section in Geneva is o"tot NEVA = 42.5 ± 1.6 nb, where the 
uncertainties are given by integrating over the different profile scale variations. The central 
value is 3.8% low and agrees within the uncertainties of ±3.8%. The uncertainty in Geneva 
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Figure 3. The Geneva partonic NLI/+LO3 result is shown compared to the analytic resummation 
of 72 matched to fixed order at NLL'+L03 in the (a) peak, (b) transition, and (c) tail regions. Also 
shown for comparison is the pure resummed result at NLL' and the fixed-order LO3 contribution. 
Figures (d), (e), and (f) show the Geneva partonic result at NNLL'+NL03 compared to the analytic 
resummation of T2 matched to fixed order at NNLL'+NLO^. The pure NNLL' resummation and 
fixed-order NLO3 result are also shown for comparison. 

that comes from integrating the spectrum over 7i > as in eq. (2.4) is much larger than the 
fixed-order uncertainty. The reason is that at any given point in the spectrum, but especially 
in the peak region, the relative uncertainties, reflecting both shape and normalization, are 
larger than in the total cross section. Hence, when integrating the spectrum to obtain the 
total cross section the uncertainties in the spectrum must cancel each other, meaning there is a 
negative correlation in the uncertainties between different regions in the spectrum. When the 
resummation and matching to fixed order is performed for the spectrum, this correlation and 
cancellation is numerically not exact for the total cross section. This is a well-known limitation 
of analytic resummation [49]. In fact, the result from Geneva is completely consistent with 
the inclusive cross section obtained using the analytic resummed result in eq. (2.4) with 
cr(T2 nt ) calculated at NNLL/+LO3 and the da/d7i calculated at NNLL'+NL0 3 . In this 
case, with T£ ut = 1 GeV, the central value is 2.9% low with uncertainties of ±3.7%. One 
way to solve this problem would be to enforce a (highly nontrivial) constraint on the profile 
scale variation to reproduce the required correlation exactly, in which case also the total 
cross section would come out exactly right. In practice, a simpler way to enforce this is to 
compute the result for the resummed cumulant as the difference between total cross section 
and integrated resummed spectrum. (This is similar in spirit to the method proposed in 
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Figure 4. The 7i distribution at NNLI/+NLO3 from Geneva before and after showering with 
Pythia 8 in the (a) peak, (b) transition, and (c) tail regions of the distribution. The analytic 
resummed result at NNLL/+NLO3 and the fixed-order NLO3 contribution are shown for comparison. 

Refs. [29, 36].) Since our focus in this paper is the differential spectrum, which serves as the 
primary input to the Monte Carlo, rather than the total cross section, we leave this for future 
improvement. 

3.2.2 Showered Results 

Next we validate our interface with the parton shower. In figure 4, we compare the NNLL'+NL03 
partonic and showered Geneva predictions with 72 CUt smeared between 0.5 — 1 GeV. We also 
show the analytic resummed NNLL/+NLO3 and pure fixed-order NLO3 spectra for compari- 
son. Before showering, the cumulant dc/d<l ) 2(72 CUt ) is in the T2 = OGeV bin and we see the 
effect of the smeared 72° ut on the spectrum in the Geneva partonic histogram in figure 4(a). 
The parton shower generates emissions inside the 2-jet bin, which fills out and determines the 
shape of the Geneva showered result in the region below 7~2 Ut and agrees remarkably well 
with analytic resummed spectrum below the cut. We can see this explicitly in figure 5 from 
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Figure 5. The peak region of the T2 distribution from Geneva after showering with Pythia 8. The 
contribution from showered events originating from 2-, 3-, and 4-parton events is shown along with 
their sum, including scale variation (blue band histogram). 

the separate contribution of 2-, 3-, and 4-parton events after showering for the central value 
in the peak region. The shape of the 2-parton showered histogram is determined by Pythia 
and the area under the histogram is the cumulant dcx/d < I ) 2(7^ cut ) calculated at NNLL/+LO3. 
The relative contribution of 3-parton and 4-parton events is determined by 73 CUt = 2 GeV, 
for which the 4-parton contribution is well behaved, giving 15% of the total cross section 
and no large cancellation with 3-parton events. These contributions all combine smoothly to 
generate the total Geneva showered result. 

The action of the shower on 3-parton and 4-parton events, which make up the spectrum 
above T™ 1 , is restricted to not change T2 by more than a power suppressed amount A7i, 
as discussed in section 3.1.3. This controls the allowed shift from the Geneva partonic to 
showered histograms in figure 4. We can see that there is excellent agreement, including 
uncertainties, between the two in the peak and transition regions. This validates that with 
our choice of A, the higher-order accuracy of the resummed 7i spectrum is not compromised 
by the shower. (Increasing A, we do observe at some point a shift of showered results away 
from partonic.) The showering does shift the T2 spectrum in the far tail away from the 
partonic result which matches the NLO3 curve, as can be see in figure 4(c). This is allowed, 
since our partonic prediction in this region becomes only leading order for 4 partons. 

3.2.3 Hadronized Results and Comparison to Data 

The full prediction for the jet resolution spectrum is obtained by turning on the hadronization 
in Pythia. This gives rise to a shift in the T2 spectrum, shown in figure 6, where "default" 
refers to the default running parameters a s {mz) = 0.1135 and Pythia e + e~ tune 1. As 
discussed in section 3.1.3, we use the standard Pythia 8 tunes without modifying any internal 
parameters. For comparison we show the Geneva hadronized result for tune 3 as well as 
for tune 1 with the world average value a s (mz) = 0.1184. We also show a comparison to 
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Figure 6. The showered NNLL/+NLO3 Geneva prediction with and without hadronization using 
the default values Pythia 8 e + e~ tune 1 and a s (mz) = 0.1135 compared to data from ALEPH [71] in 
the (a) peak, (b) transition, and (c) tail regions and to OPAL [72] in the peak and transition regions. 
The ratio of Geneva predictions to the ALEPH data is shown in (d). Also shown is the Geneva 
prediction at the central scale with a s (mz) — 0.1184 and e + e~ tune 3. 



experimental data from ALEPH [71] and OPAL [72]. We only show ALEPH data in the tail, 
since the OPAL data in this region is sparse. These measurements are fully corrected to the 
particle level, allowing us to directly compare to our hadronized predictions. Since the data 
are normalized to the total cross section, we rescale them to the total NNLO cross section 
and convert from thrust T to 7i = E cm (l — T). This allows us to directly compare the data 
to the absolute cross section predictions in Geneva, unlike a comparison between normalized 
spectra which would only test the shape. The Geneva prediction at the default values agrees 
impressively well with the data within uncertainties across the peak and transition regions 
and into the tail. The difference in the far tail is expected since here fixed-order contributions 
beyond LO4 are important, which are not yet included in our results. 

The partonic Geneva prediction does not include nonperturbative effects in the soft 
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function of 0(Aqcd/72), nor power corrections of the form 0(AQcn/E cm ). Since we strongly 
constrain the action of the Pythia parton shower to not change the analytic resummed 
NNLL/+NLO3 result as discussed in section 3.1.3 and demonstrated in figure 4, we expect the 
hadronization in Pythia to supply these missing nonperturbative effects. In effect, Pythia 
provides a well-tested model of the nonperturbative soft function in eq. (3.9). We show the 
hadronized Geneva result with Pythia e + e~ tune 3 at the central scale in figure 6 as a 
measure of the uncertainty from hadronization. Tune 3 turns out to give a smaller shift due 
to hadronization than tune 1, which makes a significant difference in the peak below < 3 GeV, 
where nonperturbative corrections are 0(1) and depend on the details of the hadronization 
model. In the transition and tail regions, we see a smaller difference with tune 3 being 
systematically lower than tune 1. This is consistent with the fact that the transition and tail 
regions are sensitive to only the first nonperturbative power correction in the soft function of 
0(A QC d/T 2 ). 

There is an important interplay between the effect of hadronization and the value of 
ots(jnz)i as discussed in Ref. [49] where a simultaneous fit to a s (mz) and the first nonper- 
turbative correction to the soft function of 0(Aqcd/72) was carried out. Generically, larger 
nonperturbative corrections shift the partonic spectrum to larger values of 75, while a smaller 
value of a s (mz) shifts the 2-jettiness spectrum downwards. This gives rise to compensating 
effects. Since tune 3 gives a smaller shift due to hadronization than tune 1, the combination 
of tune 3 with a s (mz) = 0.1135 gives an estimate of the lower bound on the combined un- 
certainty of these two effects in the transition and tail region, while the combination of tune 
1 and a s (mz) = 0.1184 gives an estimate of the upper bound. This is illustrated very well 
in the ratio of Geneva to ALEPH data in figure 6(d). Both are, however, still within the 
perturbative uncertainties from Geneva across most of the transition and tail regions. 

We have also checked that the nonperturbative shift from Pythia tune 1 is of similar 
size as expected from the fit results in Ref. [49]. This is consistent with the fact that it gives 
a good description of the data when used together with their fitted value of a s {mz)- Hence, 
we use tune 1 with a s (mz) = 0.1135 as the default since it agrees best with the data in the 
peak and provides a consistent description of the data across larger values of T2 ■ 

3.3 Predictions for Other Event Shapes 

In this section we present Geneva's predictions for a variety of dijet event shape variables. 
Examining observables other than the jet resolution variable we use as input serves to validate 
our master formula at the fully differential $2,3 level [see eq. (3.2)] rather than its projection 
onto the Ti spectrum [see eq. (2.24)]. Event shapes are particularly useful to consider because 
there exist both higher-order resummed results and precision LEP data with which we can 
compare. 

By construction Geneva correctly predicts other observables at NLO3, while maintaining 
the correct inclusive cross section. However, as discussed in section 2.2.3, it is an important 
open question to what extent the NNLL^— I-NLO3 resummation of the T2 spectrum increases 
the accuracy of resummed predictions for the other observables (beyond the partial NLL order 
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naively expected by interfacing with the parton shower) . While other observables will not be 
predicted at the same resummed order as T2 , the accuracy of the predictions for event shapes 
is expected to increase as a function of their correlation with 2-jettiness. The comparison 
of Geneva to the higher-order analytic resummation of event shapes plays a crucial role in 
numerically testing the accuracy achieved in our approach and validating the event-by-event 
perturbative uncertainties we predict. 

We present results for the C-parameter [73-75], heavy jet mass (p) [76, 77] and jet 
broadening (B) [78, 79] event shapes. These are defined as follows: 



^2 pk ) for * = 1 > 2 > 

(3.20) 

where fix is the thrust axis and is used in heavy jet mass to divide the event into two 
hemispheres, hemii^, with respect to which the masses M\p are measured. C,p, and B 
provide a useful range of event shapes to compare to since their resummation structure is 
increasingly different from that of 7i- The resummation of C-parameter is precisely the 
same as 7i to NLL and has the same convolution structure as eq. (3.8) beyond. Heavy jet 
mass has a different convolution structure from 7i- Both p and 7i are projections of the 
same doubly differential spectrum dcr/dM^dM|, where 7i is related to the sum and heavy 
jet mass to the maximum of the hemisphere masses. Of the event shapes we consider, jet 
broadening is most different from T2] it measures momentum transverse to the thrust axis 
and in the dijet limit is sensitive to the recoil of the thrust axis due to soft emissions [80], 
unlike 7i- This complicates the higher-order resummation of jet broadening, which was only 
recently extended to NNLL# [81] and gives a logarithmic structure that is very different to 
7i- As a result, jet broadening provides a highly nontrivial test of the accuracy and theory 
uncertainties of the Geneva prediction. 

For each of these observables we compare to analytic resummed predictions as well as 
the NLO3 fixed-order contribution from Event2. We present new results for the analytic 
resummation of C-parameter at NNLL^+NL03, extending the previous NLLc resummation 
of [82] 12 . Note the subscript on the order of resummation indicates the observable for which 
analytic resummation was carried out. Since resummed results for jet broadening do not 
exist at NNLL' B , we compare to the highest available resummation NNLLB+LO3, where we 
use the results of [81], which we extend to include fixed-order matching that is necessary to 
describe the tail and transition regions. Finally for heavy jet mass N 3 LL p resummed results 

12 We thank Vicent Mateu and Iain Stewart for pointing out to us the relationship between thrust and 
C-parameter in SCET. 
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exist [48], however we show the NNLL^+NLOs resummation since this is consistent with the 
highest I2 resummation we use. 

It is important to note that all running parameters were set based on the T2 spectrum 
alone and no further optimization was carried out for other observables. This ensures that 
our results for other observables are true predictions of the Geneva framework. 

3.3.1 C-parameter 

In figure 7 we show the Geneva prediction for C-parameter both at the partonic level and 
showered, using NLL^+L03 resummation as input to our master formula in figures 7(a) and 
7(b) and at next higher order NNLL^- I-NLO3 in figures 7(c) and 7(d). We compare this to 
the analytic resummed C-parameter prediction at the same order as the 7i resummation we 
input, as well as one order lower. The comparison of the Geneva prediction at different orders 
in the peak and transition regions is useful because it highlights the features of resummation 
that are consistently captured by our implementation. In the tail region, figure 7(e), where 
the comparison to the NLO3 fixed-order result is most relevant, we only show our highest 
order NNLL^+NL0 3 Geneva result. 

We see the effect of the cut on 2-jettiness of T 2 cnt = 0.5 — 1 GeV up to C = 0.066 in the 
partonic prediction from Geneva in figures 7(a) and 7(c) since C < &T2/Q [82]. Interfacing 
with the shower generates emissions inside the jets and fills out the region below C = 0.066. 
The action of the parton shower is restricted based on the constraints on 7i discussed in 
section 3.1.3. This effectively constrains the C-parameter distribution as well, giving very 
little change from the partonic to showered predictions at both resummation orders, except 
in the multijet region of the far tail. Here the constraints on the shower are looser, reflecting 
the fact that our prediction is correct at LO4. The size of the shift from the partonic to the 
showered result in the peak and transition regions is a measure of the correlation between the 
C and 7i event shapes, where although the two differ beyond NLL, their logarithmic structure 
is the same. It is worth noting however, that despite the similarity of the resummation 
structure between C and 7i, the shape of the C-parameter spectrum is very different, with 
the singular terms dominating the nonsingular for a much larger region of the spectrum. 

One might naively expect that the accuracy of the resummation achieved in Geneva for 
any observable other than 7i would only be the partial NLL of the parton shower. However it 
is clear that the Geneva prediction at NLL^+L03 in the peak and into the transition region, 
figures 7(a) and 7(b), is much more consistent with NLL'.+LOs than NLLc resummation, 
both in its central value and the size of the perturbative uncertainties it predicts. This appears 
to hold even in the peak region below C ~ 0.05 where the parton shower determines the shape 
of the spectrum. Going to one higher order in figures 7(c) and 7(d) we see that the same 
pattern holds: the Geneva prediction is consistent with the higher-order NNLL^.+NL03 
resummation rather than NLL^.+L03, including uncertainties. This is particularly clear 
in the peak region where the central values of the two analytic resummation orders are 
significantly different and Geneva tracks the NNLL^,+NL03 prediction. The convergence of 
the Geneva result for C-parameter from NLL^+L03 to NNLL^+NL03 demonstrates the 
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consistency of the Geneva implementation including the event-by-event uncertainties for this 
observable. Although the accuracy of the Geneva prediction for C-parameter is not formally 
of the same order as the T2 resummation we used as input to the master formula, the fact 
that it matches the analytic C-parameter resummation remarkably well both at NLL^+LOs 
and at NNLL^+NLOs shows that numerically the accuracy achieved is very close. 

The Geneva uncertainties in the transition region start to shrink relative to the analytic 
resummation as we interpolate to the fixed-order NLO3 result. In the tail region the par- 
tonic Geneva prediction matches smoothly to the fixed-order NLO3 result past the Sudakov 
shoulder at C = 0.75 [83], demonstrating the validity of the multiplicative implementation of 
da"/d$3(72) in eq. (3.3) in this limit. 

The Geneva prediction including hadronization with the default running values of Pythia 
e + e~ tune 1 and a s (mz) = 0.1135 is shown in figure 8 compared to ALEPH and OPAL data 
rescaled to the NNLO inclusive cross section. Geneva agrees with the data remarkably 
well across the entire distribution up to the multijet region in the tail. We show the ef- 
fect of a s (mz) = 0.1135 with tune 3, which gives a smaller correction from hadronization 
than tune 1, as seen from the size of the shift from the Geneva unhadronized result to the 
hadronized in figure 8. We also show the Geneva prediction at the central scale using the 
world average a s (mz) = 0.1184 and tune 1. These two combinations provide an estimate of 
the upper and lower bounds on the combined uncertainty of the nonperturbative effect and 
ot-sijnz) value in the transition and tail regions, as discussed in section 3.2. The ratio of the 
Monte Carlo to data in figure 8(d) shows that they are both largely within the perturbative 
uncertainties from Geneva in these regions. Of the values we consider, the default tune 1 
with a s (mz) = 0.1135 gives the best agreement with the data across the C spectrum and is 
consistent with our findings for the 7i distribution. 

3.3.2 Heavy Jet Mass 

The Geneva prediction for heavy jet mass is shown in figure 9, where we compare the partonic 
and showered results using NLL^— I-LO3 resummation in the master formula in figures 9(a) 
and 9(b) to the analytic resummation of p at NLLp and NLL^-(-L03. In figures 9(c) and 9(d) 
we show the same results at one order higher, comparing NNLL^+NL03 Geneva results to 
NLLp+L03 and NNLL^+NL03 analytic p resummation. In the tail we show only the highest 
order Geneva and resummed results along with the pure fixed-order NLO3 contribution, 
since this is sufficient to demonstrate the behavior in this region. 

In the peak region, figures 9(a) and 9(c), we see the effect of 7~2 Ut on the partonic p 
spectrum up to p = 1 GeV/E cm = 0.011, which is smoothly filled out by interfacing with the 
parton shower. The Geneva showered prediction in figure 9(a) shows impressive agreement 
with the NLLp+L03 resummed result in the peak region, including uncertainties. The im- 
provement in accuracy of the Geneva prediction for heavy jet mass over the partial NLL 
resummation provided by the parton shower is clear. Going to one higher order in figure 9(c) 
the Geneva prediction is more consistent with the NNLL^+NLOs, with which it agrees 
within uncertainties, rather than NLLp+L03 result. The perturbative uncertainties of the 
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Geneva showered prediction are larger than those at NNLL^+NL03 and smaller than at 
NLL^+L03. This is consistent with the fact that while the Geneva prediction for p is not 
formally of the same order as the 7i resummation that is input into the master formula, there 
is a gain in accuracy going from Geneva at NLL^+L03 to NNLL^+NLC>3 that is transferred 
to the p prediction. 

In the transition region, both at NLL^- I-LO3 and NNLU7-+NLO3, adding the parton 
shower gives rise to a larger shift from the partonic spectrum than for the C-parameter, 
because heavy jet mass is less correlated with T2 than C. This shift is necessary to obtain 
agreement with the NNLL^+NL03 resummation within uncertainties in figure 9(d). The 
partonic Geneva prediction in this region is more consistent with NLL^+L03 analytic re- 
summation, which is higher than the NNLL^+NLOs result. By restricting the change in 7i 
due to the shower, we are constraining the sum of the hemisphere masses, in eq. (3.20). 
For a given value of 7i, p is largest when either hemisphere mass is zero and so p = 7i- 
Adding the parton shower tends to make this mass nonzero (while constraining the sum) and 
therefore gives an overall shift of the spectrum to lower values of p. In the tail region the 
partonic spectrum interpolates to the fixed-order NLO3 result, as expected, with the shower 
giving rise to a larger shift in the multijet region where our constraints are looser. 

In figure 10 we show the showered Geneva prediction with and without hadronization, 
with our default parameters. As before we compare to ALEPH and OPAL data, which shows 
impressive agreement with the data within uncertainties across all 3 regions of the p spectrum, 
with the expected deviation in the multijet region of the far tail. As discussed previously, 
tune 3 with a s (mz) = 0.1135 and tune 1 with a s (mz) = 0.1184 provide bounds on the 
estimate of the combined uncertainty from these two inputs. It is interesting to note that 
heavy jet mass is relatively insensitive to hadronization in the transition and tail regions. This 
is demonstrated by the shift from the shower only to the full hadronized result in figures 10(b) 
and 10(c), as well as the small change in the default central value from using tune 1 to tune 
3 above p ~ 0.1 seen in figure 10(d). This breaks the coupling in some respect between the 
nonperturbative effects and the value of a s in this region and suggests that a s (mz) = 0.1135 
provides better agreement with the data. 

3.3.3 Jet Broadening 

Finally we turn to the results of Geneva for jet broadening, which is the most orthogonal 
event shape to our jet resolution variable that we consider. In figure 11 we show the Geneva 
partonic and showered results using NLL^— I-LO3 resummation in figures 11(a) and 11(b) and 
NNLL^+NL03 resummation in figures 11(c) and 11(d). We compare these to the analytic 
NLL# and the best available NNLLB+LO3 resummed prediction. Note that we would like 
to compare the NLL^+L03 resummation in Geneva to the resummation of B at the same 
order. However, since going from NLL^ to NNLL^, which incorporates a 2 s \n.B terms into 
the resummation, is a comparatively small effect in this case, we will find it useful to compare 
the NLL^— I-LO3 prediction with the analytic NNLLB+LO3 result. In the tail we compare 
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to the fixed-order NLO3 result from Event2, which is the more relevant comparison in this 
region. 

The effect of the cut on T2 extends up to B ~ 0.055 in the peak region and is smoothly 
removed by attaching the parton shower. Both at NLL^+L03 and NNLL^- I-NLO3 there is 
a significant shift induced by the parton shower across the jet broadening spectrum towards 
larger values of B. The size of this shift is a measure of the the lack of correlation between B 
and 72, the variable used to constrain the parton shower, and is therefore progressively larger 
for C, p, and B as we have seen. As discussed in section 2.2.3, in the IR limit where both 
T2, B — > one might expect to see some improved accuracy of the Geneva prediction over the 
partial NLL of the parton shower, since the higher-order resummation of 72 provides a better 
description in this region. This is consistent with figure 11(a) where the showered Geneva 
prediction agrees well with the NNLL5+LO3 resummed result, including uncertainties in the 
peak region. In the transition region in figure 11(b), the central value of the Geneva showered 
prediction agrees with the NNLL5+LO3 resummed result within uncertainties; however, the 
uncertainties from Geneva are smaller than the corresponding analytic ones, which suggests 
that in this region they may be underestimated. 

In the far transition region and into the tail, the uncertainties in Geneva generically are 
smaller than the corresponding analytic resummation and of order the NLO3 scale variation, 
as seen for example in the Ti, C, and p spectra. This difference arises because Geneva 
multiplicatively interpolates to the fixed-order result, while the analytic resummation does 
an additive matching, as discussed in section 3.2.1. For an observable such as jet broadening, 
the lack of correlation with 7i means that larger values of the 72 spectrum contribute at 
smaller values of B. This can lead to an underestimate of the uncertainties from Geneva at 
intermediate values of B where the resummation is still important. 

Going to higher order in Geneva in figures 11(c) and 11(d) the uncertainties of the 
showered prediction decrease and overlap with the NNLLS+LO3 uncertainties over much of 
the peak and transition regions. It would be interesting to compare the Geneva prediction 
to the next higher order analytic resummed jet broadening prediction to numerically test 
the accuracy achieved, however this is not yet available. Determining the formal accuracy of 
the Geneva prediction for a given observable and systematically including the uncertainty 
associated to the lack of correlation with 72 are next steps that we leave for future work. 

The partonic Geneva jet broadening prediction interpolates smoothly to the fixed-order 
NLO3 result in the tail, with uncertainties that match the fixed-order result in this region. 
As for other observables, this validates the behavior of the fully differential master formula 
in eq. (3.3) in this limit. 

In figure 12 we show the hadronized Geneva prediction for jet broadening compared 
to data from ALEPH and OPAL, which shows good agreement within uncertainties across 
the peak and transition regions, and is low as expected in the far tail. The uncertainty 
from the Pythia tune and value of a s are indicated by the central values of the tune 3 
and a s = 0.1184 histograms, which agree within the perturbative uncertainties of the default 
Geneva prediction across most of the transition and tail regions. As for the other observables 
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we see better agreement in the peak (below B ~ 0.1) with Pythia tune 1 and a s = 0.1135. 
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Figure 7. The C-parameter partonic and showered Geneva predictions are shown compared to 
the analytic resummation of C-parameter at different orders. The Geneva result at NLL?r-+L0 3 is 
compared to NLLc and NLLp-|-L0 3 in (a) and (b). In (c) and (d) the Geneva prediction at one 
order higher, NNLL/ r +NL0 3 , is compared to NLL^+L0 3 and NNLL^,+NL0 3 , while in the tail (e) 
we also show the fixed-order NL0 3 prediction from Event2. 
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(c) Tail Region (d) Ratio of Geneva to Data 



Figure 8. The C-parameter distribution comparing Geneva with and without hadronization using 
Pythia 8 e + e~ tune 1 and a s (mz) = 0.1135 is shown compared to ALEPH data in the (a) peak, (b) 
transition, and (c) tail regions and to OPAL data in the peak and transition regions. The ratio of the 
Geneva predictions to ALEPH data is shown in (d). Also shown are the Geneva predictions at the 
central scale with a s (mz) = 0.1184 and e + e~ tune 3. 
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Figure 9. The heavy jet mass partonic and showered Geneva predictions are shown compared to 
the analytic resummation of p at different orders. The Geneva result at NLL^+L0 3 is compared 
to NLLp and NLL^+L0 3 in (a) and (b). In (c) and (d) the Geneva prediction at one order higher, 
NNLL^+NL0 3 , is compared to NLL/ +L0 3 and NNLL/,+NL0 3 , while in the tail (e) we also show 
the fixed-order NL0 3 prediction from Event2. 
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Figure 10. The heavy jet mass distribution comparing Geneva with and without hadronization using 
Pythia 8 e + e~ tune 1 and a s (mz) = 0.1135 is shown compared to ALEPH data in the (a) peak, (b) 
transition, and (c) tail regions and to OPAL data in the peak and transition regions. The ratio of the 
Geneva predictions to ALEPH data is shown in (d). Also shown are the Geneva predictions at the 
central scale with a s (mz) — 0.1184 and e + e~ tune 3. 
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Figure 11. The jet broadening partonic and showered Geneva predictions are shown compared to 
the analytic resummation of B at NLL B and NNLL B +L0 3 . The Geneva result at NLL^+L0 3 is 
shown in (a) and (b), and at one order higher, NNLL^+NLOs, in (c) and (d). In the tail (e) we also 
show the fixed-order NLO3 prediction from Event2. 
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(a) Peak Region (b) Transition Region 




Figure 12. The jet broadening distribution comparing Geneva with and without hadronization using 
Pythia 8 e + e~ tune 1 and a s (mz) = 0.1135 is shown compared to ALEPH data in the (a) peak, (b) 
transition, and (c) tail regions and OPAL data in the peak and transition regions. The ratio of the 
Geneva predictions to ALEPH data is shown in (d). Also shown are the Geneva predictions at the 
central scale with a s (mz) = 0.1184 and e + e _ tune 3. 
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4 Application to Hadronic Collisions 

In this section, we apply the framework developed in section 2 to hadronic collisions and 
present first results from the implementation in the Geneva Monte Carlo. To accommodate 
hadrons in the initial state, special consideration is required for each component of our master 
formula. Our goal is to demonstrate that the same methods can be applied in a hadronic 
environment to obtain a consistent description at the next higher perturbative accuracy. We 
use Drell-Yan production pp — > Z/j — > £ + l~ as a concrete example to study the framework, 
deferring a detailed comparison with Tevatron and LHC data to later work. 

In hadronic collisions iV-jettiness can be used as a jet resolution variable, with the ob- 
servable taking the initial states into account. The theoretical framework exists to resum 
iV-jettiness at hadron colliders, and this resummation has been applied to several processes 
[66, 84, 85]. Similarly, the techniques required to perform the next-to-leading order calcula- 
tions at hadron colliders are known. A phenomenological study additionally requires Geneva 
to be interfaced with a parton shower and hadronization model that includes multiple parton 
interactions. 

In the Drell-Yan example the 0-jet resolution variable is beam thrust, 7o, which is the 
analog to 2-jettiness, 7i, used in the e + e~ application 13 . The resummation of beam thrust is 
performed to NNLL, and the jet multiplicities at fixed order are calculated at NLOo and LOi. 
The prediction of Geneva is compared to the analytic resummation of To at NNLL+LCq 
as a demonstration that the program correctly describes the matching between 0- and 1-jet 
multiplicities. Finally, we discuss how the accuracy of these ingredients can be improved and 
the challenges present in applying Geneva to hadron collisions. 

Beam thrust is defined as a sum of contributions from particles in the final state [43]: 

% = ^mm{n a -p k ,n b -p k }, (4.1) 

k 

where the observable is evaluated in the center-of-mass frame of the hard partonic collision. 
The n a fi are light cone vectors along the beam (I) axis, with n a = (1,-2) and n b = (lj — z). 
Beam thrust can be evaluated in any frame by performing a longitudinal boost on eq. (4.1). 

With the addition of more final state jets, the iV-jettiness definition can be generalized 
from the 0-jet case: 

Tn = ^min{n a ■ p k , n b ■ p k , m ■ p k , . . . , n N ■ p k ) . (4.2) 

k 

As for beam thrust, this observable is evaluated in the partonic center-of-mass frame. The 
m = (l,Wi) for i = 1, ...,7V are light cone vectors along the jet directions. Note that the 
above definition of T/v is a simple extension of the observable for e + e~ collisions, which has 
no contribution from the beam directions but is otherwise identical. 

13 Since we are mainly interested in QCD corrections, we have chosen the subscript on T to indicate the 
multiplicity of jets in the final state. 
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4.1 Master Formula and Ingredients for Hadronic Collisions 

As in the e + e~ case, the master formula for the cross section in Geneva is given by eq. (2.27), 
eq. (2.29), and eq. (2.36). To match the 0- and 1-jet multiplicities for a general process, the 
master formula is 

where the 0-jet cumulant, d<7/d<I>o(7^| CUt ) and the 1-jet spectrum, da/d&i(To), are 

W° da 1 ' 



j J rcsum 

■Uo) 

1 



d$ ' d$ 



FO. 



d$i d$ 



d(J rosum / da 



rcsum 



d$ d75/ d$ dT 



FO. 



(4.4) 



$0 is the phase space for the hard scattering that produces the 0-jet final state, and <E>i 
includes the additional phase space for the final state jet. 

In Geneva, the contributions to these cross sections are calculated separately for each 
parton subprocess. Because the fixed-order matching in the 0-jet cumulant is performed 
additively, the net weight in the 0-jet cumulant after summing over events is the same as 
the flavor-summed cumulant. In the 1-jet spectrum the fixed-order matching is performed 
multiplicatively, meaning the sum over all events for a given 7o has a different cross section 
than if we used flavor-summed components for the different terms in the matching formula. 
The two approaches agree up to higher order corrections, but the former approach is natural 
in the Monte Carlo. In the following subsections, we will discuss how the resummed and 
fixed-order contributions to the master formula are obtained. 

4.1.1 Resummation 

Like 7i for e + e~ collisions, beam thrust can be factorized in SCET and the resummation 
can be carried out systematically to higher orders. The factorized beam thrust spectrum for 
Drell-Yan is given by [31] 



da das 



^Y^Hij(Q,ix) j dt a dt b B i (t a ,x a ,ii)Bj(t b ,x b ,ii)Si j (To - t ^^- ,fxj , (4.5) 



d$ dT d$ 

where Q is the dilepton invariant mass, and $0 is the phase space for the qq —> £ + £~ hard 
scattering. The momentum fractions x a ^ are defined in terms of the total rapidity Y of the 
final state from the hard scattering, 

Q Y _ Q -Y 



„ ,, e Y , x b = ^e~ Y . (4.6) 

"cm "cm 

Comparing eq. (4.5) to the e + e" analog, eq. (3.8), it is clear that the factorization theorems are 
structurally similar. The chief difference is that while the jet functions in eq. (3.8) parametrize 



- 51 - 



the collinear evolution of final state jets, the beam functions in eq. (4.5) parametrize the 
collinear evolution of the incoming partons as well as the nonperturbative process of extracting 
high energy partons from the proton. The beam functions can be further factorized into 
a convolution between the parton distribution functions, fj, and perturbatively calculable 
Wilson coefficients, Xij [31, 86] 



Note that due to initial state radiation, the x a ,b are distinct from the Bjorken variable £ 
appearing in the convolution above that give the momentum fraction of the energetic partons 
that are liberated from the proton. The sum over partons i,j in eq. (4.5) is a sum over 
flavor singlet quark- antiquark combinations, such as uu or bb. For each flavor, the beam 
functions are different, and the hard function Hij differs for up- and down- type quarks due 
to the different electroweak couplings with the gauge boson. This flavor sum is an important 
consideration when implementing the master formula in Geneva, since the Monte Carlo 
generates events for each flavor combination, and the flavor sum is performed in the sum over 
events. 

For processes with final state jets, the extension of the beam thrust factorization theorem 
to iV-jettiness is known and has the schematic form [43, 44] 



The trace is over the non-trivial color structures that can exist in the hard and soft functions. 
Additionally, there is a sum over parton channels for the hard scattering that are labeled by 
the index k. The additional jets are associated with additional collinear sectors in SCET, and 
the factorization theorem reflects this by including additional jet functions. The soft function 
also changes to account for the soft radiation between the final-state jets and the initial-state 
radiation from the colliding partons. 

The factorization theorems in eqs. (4.5) and (4.8) can be used to perform the resummation 
for both the spectrum and cumulant. Although these factorization theorems directly describe 
the spectrum in 7o or 7~n, they can be integrated over the observable to obtain the cumulant. 
The perturbative part of each function in the factorization theorem is calculable, and for 
many processes the functions are known to high order. Each function is associated with a 
scale that is connected to the physical degrees of freedom that the function describes. As in 
the e + e~ case, renormalization group evolution resums the large logarithms of ratios of these 
scales (see table 2). 

4.1.2 Fixed Order 

Following eq. (4.3) we need to define the 0-jet cumulant do"/d<I>o(7^ cut ) and the 1-jet spectrum 
dcr/d<3?i(7o)- At the order we are interested in, the 1-jet spectrum will be given by the tree- 
level cross section B\(&i) for the process pp —> Z/~f* — > l + £~j. The 0-jet cumulant is given 




(4.7) 





by 
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where the corresponding parton distribution functions have been included into the definitions 
of the Born, virtual, and real emission cross sections, 

B N ($ N ) = f a (x a , n F )f b (x b , ij, f )B n (® n ) , 

V n ($n) = fa{x a , fiF)fb{xb, Hf)V n {$ n ) . (4.10) 

In addition, in order to account for the incomplete cancellations of initial-state collinear 
singularities, we have included one collinear counterterm G a ^ for each initial state parton, 

G a ($ N ) = fa(x a , fl F )f b (x b , HF)Ga{$N) , 

G b ($ N ) = f a (x a , fi F )f b (x b , HF)Gb{<S>N) ■ (4.11) 

Assuming the UV divergences of Vo have already been taken care of by a proper renormal- 
ization procedure, the remaining divergences present in B\, Vo, and G a ,b are of IR origin. We 
handle these divergences with the FKS subtraction procedure. After having partitioned the 
phase space into nonoverlapping regions m which contain at most one collinear and one soft 
singularity, by means of a set of ^-functions, the final formula, including the subtraction 
counterterms, is 



da 
d$n 



cut\ 



Bo(*o) + Vo('I>o) +/($(,) + 
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d$ m 
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rad 



rad 
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rad 



rad 



d$! 



rad 



+ J d$ rad , a G (*i, B ) + j d$ radi & G 6 (*l,6) • 



(4.12) 



As mentioned in section 3.1.2 we choose to partition the phase space by means of 
^-functions that depend on the jet resolution variable. It is therefore crucial to evalu- 
ate the jet resolution variable and the subtraction in the same frame. This ensures the proper 
cancellation of IR singularities by subtraction counterterms. The preferred frame for the 
fixed-order calculations is the partonic center-of-mass frame, since the subtraction is most 
naturally expressed in terms of variables defined in that frame. Also, the jet resolution vari- 
able, eq. (4.1), is defined in this frame and resummation can be performed in it. Therefore, 
our approach is to perform the entire calculation in the partonic center-of-mass frame. 

At this point all the ingredients of eq. (4.12) are known and available in the literature [17, 
67]. Note that additional complications arise when extending this construction to higher 



- 53 - 



multiplicities, because one must use a map that preserves the value of To (up to power 
corrections), just as for T2 in the e + e~ — > 3 jet case. For our pp — > Z/7* — > £ + £~ study, 
this problem can be avoided, since we currently include only up to one extra parton and, 
consequently, the 7o(^i) value is well-defined. In order to obtain To(&n) with N > 1 in 
general, this issue may be addressed in a similar fashion to what has been done for e + e _ . 

4.2 Application to Drell-Yan Production 

We study Drell-Yan production in pp collisions at E cm = 8 TeV in Geneva, sampling the 
invariant mass Q of the £ + £~ pair around the Z pole, between Mz — WTz and Mz + lCT^, 
where M z = 91.1876 GeV is the mass of the Z and its width is T z = 2.4952 GeV [54]. The 
dominant contribution in this range of Q comes from Z exchange, although the photon does 
contribute. Profile scales identical to those used in the e + e~ 7i resummation are used, which 
is justified since the logarithmic structure is the same between the two observables. The 
resummation is turned off just above To ~ M^/2, and for greater To the spectrum reproduces 
the fixed-order distribution. 

In figure 13, we show the analytic beam thrust resummation at NLL and NNLL+LOi 
in the peak, transition, and tail regions. In the peak and transition regions the resummed 
result converges well. At the end of the transition region and in the tail region, the pure 
NLL resummed distribution goes to as the resummation is turned off but the NNLL+LOi 
distribution moves into the fixed order result.. 

Implementing the Drell-Yan process in Geneva allows us to study the feasibility of the 
multiplicative matching for the spectrum in eq. (4.4) and compare with the analytic resummed 
distribution. We show this comparison in figure 14, where the analytic curve is evaluated at 
NNLL+LOi. Additionally, we show the NNLL and LOi distributions separately. Overall, the 
partonic Geneva distribution is quite close to the NNLL+LOi distribution, both in terms of 
the central values and the size of uncertainties. In the peak region, the spectra match closely 
and agree well with the pure NNLL resummed distribution. At the low end of the transition 
region, the resummed spectra are still in fair agreement while moving to higher To values 
the Geneva partonic prediction and the NNLL+LOi distributions begin to systematically 
deviate from the NNLL distribution. This deviation arises from the LOi nonsingular terms 
that are present in the matched spectrum but absent in the pure resummed one. In the tail 
region, the Geneva partonic and NNLL+LOi predictions move into the LOi spectrum. After 
the resummation has been turned off, these spectra match the LOi precisely in central value 
and uncertainties, as expected. Note that in the corresponding comparison in the e + e~ case, 
shown in figure 3(a)-3(c), the analytic and the partonic Geneva distributions are in closer 
agreement because the resummed components of the multiplicative matching in eq. (2.36) 
include the nonsingular terms at LO3, which are known analytically. These are not included 
in the Drell-Yan case, and so the difference between the analytic and Geneva distributions in 
figure 14 is more sensitive to the subleading corrections that the nonsingular terms generate. 

As in the e + e~ case, the uncertainty bands for the resummed curves and Geneva predic- 
tions in figures 13 and 14 are obtained by adding the fixed-order and resummation uncertain- 
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Figure 13. Analytic resummation of 7o matched to fixed order in the (a) peak, (b) transition, and 
(c) tail regions. The central value is shown along with the band from scale uncertainties, as discussed 
in section 4.2, at NLL and NNLL'+LOi. 



ties in quadrature. In the peak and transition regions of the distribution, the resummation 
uncertainties dominate, while the fixed-order uncertainty dominates as the resummation is 
being turned off in the tail region. Comparing the uncertainty of the resummed distributions 
with that of the LOi distribution, which is much smaller, one can see that the fixed-order un- 
certainty is an underestimate of the missing higher-order terms. The reason for this is twofold: 
the missing large logarithmic corrections at higher orders, whose associated uncertainties are 
instead included in the resummed distribution, and the partial cancellation between the renor- 
malization and factorization scale dependence, whose variation are correlated in the results 
we show. 

In both the e + e~ and Drell-Yan processes, the partonic Geneva spectrum is determined 
by eq. (2.46), which for an event multiplies the fully exclusive fixed-order cross section by 
the ratio of the resummed cross section for the jet resolution variable divided by the fixed- 



- 55 - 





Figure 14. The Geneva partonic NNLL+LOi result is shown compared to the analytic resummation 
of 7o matched to fixed order at NNLL+LOi in the (a) peak, (b) transition, and (c) tail regions. Also 
shown for comparison is the pure resummed result at NNLL and the fixed-order LOi contribution. 

order expansion of that resummation. Compared to the e + e~ case, where each subprocesss 
contributing to the cross section is trivially proportional, in Drell-Yan the convolution with 
the PDFs requires treating every possible qq initial state separately, both in the fixed-order 
and the resummed cross sections. In Geneva, the flavor sum is performed in the Monte Carlo 
sense, since every event has a definite flavor for the initial-state quarks and the correct flavor- 
summed cross section is obtained after a sum over all events. This means that the separate 
factors in eq. (2.46) are evaluated for an individual flavor, and the entire expression is summed 
over flavors. In the analytic resummation, since the matching between the resummed and 
fixed-order cross sections is additive, there is instead only one way to perform the sum over 
flavors. 

A version of the master formula where both the resummed and the resummed-expanded 
are separately flavor summed before entering eq. (4.4) would be equally valid. We have 
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checked, however that this is has a very minor effect and is not the main contribution to the 
apparent differences of the Geneva predictions with the analytic resummed cross section. 
In fact, the reason of the discrepancy is in the difference of higher order terms which are 
included in the Geneva multiplicative approach with respect to the additive matching used 
in the analytical calculation. This can also be seen as an indication of the relative freedom in 
implementing the master formula in equation eq. (2.46) at a given perturbative accuracy. As 
one can evince from figure 3(d)-3(f), should the resummation and fixed-order calculations be 
evaluated at the next order in perturbative accuracy, the size of the yet-missing terms would 
decrease and, consequently, the difference between Geneva and the analytic results would 
be reduced. 

The Geneva implementation in this example can be extended to higher accuracy both 
in terms of fixed-order matrix elements and resummation. An equivalent accuracy to the 
e + e~ results shown in section 3 can be achieved if the fixed-order matrix elements for the 
1-jet multiplicity are calculated at NLO, the 2-jet multiplicity are calculated at LO, and the 
resummation of beam thrust is continued to NNLL'. Although this is beyond the scope of this 
work, we nonetheless demonstrate that the Geneva framework is capable of merging matrix 
elements of different jet multiplicities beyond the lowest order. As in the e + e~ case, jet 
multiplicities are defined using a physical jet resolution variable, which allows for a consistent 
extension of the entire framework to 0(a s ) perturbative accuracy. 

In section 3, we found that NNLL' resummation of the jet resolution variable T2, when 
combined with the parton shower, was capable of describing the spectrum in other two-jet 
observables with an accuracy that clearly exceeded NLL, the naive accuracy of the parton 
shower. The resummation of 7i captures an important set of logarithms that are correlated 
with other two-jet observables, and when combined with the fully exclusive, all-orders de- 
scription of the parton shower, the accuracy of other observables can be improved beyond 
NLL. At a hadron collider, the effective dynamic range of observables is much larger, meaning 
the correlation between the jet resolution variable and another observable of interest may be 
small. In this case, the parton shower may play a greater role in determining the spectrum, 
and hence the accuracy, of other observables. We will investigate these features in a future 
work. 

5 Conclusions 

In this paper we have shown how to combine higher-order resummation of a jet resolution 
variable with fully differential next-to-leading order calculations to extend the perturbative 
accuracy of cross sections beyond the lowest order for all values of the jet resolution scale. 
This framework has been interfaced with a parton shower and hadronization to give the 
Geneva Monte Carlo program. 

Our framework provides both the versatility of fixed-order calculations and the accuracy 
of higher-order analytic resummation. From the point of view of Monte Carlo generators, the 
Geneva approach allows the combination of higher-order resummations with higher fixed- 
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order calculations. From the point of view of resummed calculations, it allows one to obtain 
a fully differential cross section that correctly resums the jet resolution variable to higher 
logarithmic accuracy. Since this construction maintains the higher perturbative accuracy 
for all values of the jet resolution scale, it naturally allows NLO calculations of different jet 
multiplicities to be combined with one another. 

The higher logarithmic resummation of the jet resolution scale allows us to use a low 
cut on the jet resolution scale, much lower than the point where fixed-order perturbation 
theory breaks down but above the nonperturbative regime. This is a major difference to 
other approaches [26-28], where the merging scale has to be chosen much larger, such that 
a s ln 2 T cut < 1. 

In this paper we have concentrated on the theoretical construction, which is valid for 
any number of jets, and for both e + e~ and hadron colliders. We have shown that one has to 
carefully choose a jet resolution variable that is resummable to higher logarithmic accuracy. 
The presence of IR divergences led us to divide the inclusive iV-jet cross section into an 
exclusive iV-jet and inclusive (N + l)-jet cross section, and we have given expressions for both 
of these that properly combine the higher logarithmic resummation with a higher fixed-order 
calculation. 

This approach has been implemented in the Geneva Monte Carlo. As a first application, 
we have presented results for e + e~ collisions. The jet resolution variable for this case was 
chosen to be 2-jettiness, which is directly related to thrust, and we combined its NNLL' 
resummation with the fully differential 3-jet rate at NLO3. Varying the profile scales and the 
renormalization scale has allowed us to obtain event-by-event uncertainties. As a final step, 
we have interfaced our perturbative result with Pythia 8, which added a parton shower and 
hadronization to our results. The parton shower adds additional radiation beyond the highest 
jet multiplicity in Geneva, and has been restricted to only fill out the jets of the exclusive jet 
cross sections at lower jet multiplicity. Since the cut on the jet resolution variable could be 
chosen to be very small, the effect of the perturbative shower (without hadronization) is rather 
small, and different tunes in Pythia do not affect the resulting distributions significantly. 
Hadronization has a significant effect, and the difference in final results due to different 
hadronization parameters is more manifest. 

We have shown that Geneva correctly reproduces the higher-order resummation of the 
thrust spectrum, even after showering, which serves as a nontrivial validation of our ap- 
proach. Using a s {mz) = 0.1135, as obtained in Ref. [49] from fits to the thrust spectrum 
using higher-order resummation, together with tune 1 of Pythia 8, we obtain an excellent 
description of ALEPH and OPAL data. The same setup was then used to predict other event 
shape variables, namely C-parameter, heavy jet mass, and jet broadening. In all cases our 
results agree remarkably well with the explicit analytical resummations, even though only 
the thrust resummation was used as an input. This comparison shows numerically that we 
achieve a higher resummation accuracy than NLL, which is what one would naively expect 
to obtain from the parton shower. This is especially remarkable for jet broadening, where 
the resummation formula has a completely different structure from the thrust resummation. 
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Comparing our results after hadronization to the data we again find excellent agreement for 
these other observables. 

Finally, we have presented first results toward a implementation for hadron colliders in 
Geneva. Choosing the Drell-Yan process at the LHC, with beam thrust as the jet resolution 
variable, we combined the resummation of beam thrust at NNLL with the leading-order 
matrix element for the emission of an extra jet. The results from Geneva agree well with 
analytical results, which shows the applicability of the framework to hadron colliders. 

As we have shown, our theoretical framework to combine higher-order resummation with 
fixed-order matrix elements and parton shower Monte Carlos is very general, and there are 
many avenues to pursue in the future. Obvious next steps for e + e _ collisions are to include 
NLO calculations for 4 jets, which would require including the logarithmic resummation for 
3-jettiness [87], as well as additional tree-level matrix elements. In addition one can consider 
the resummation for other jet resolution variables. For hadronic collisions, next steps are to 
include the resummation and NLO calculations for higher jet multiplicities, as well as adding 
parton showering and hadronization using the different available programs. 
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